use api of vcfpp 0.3.7

Author: Zilong-Li

cran-comments.md

@@ -1,4 +1,2 @@
-1. address clang-UBSAN issue by including the latest htslib-1.19.1 (https://github.com/samtools/htslib/releases/tag/1.19)
-2. reduce size of package
-3. add copyrights and authors of htslib
-4. new function `vcfcomp`
+1. address issues on https://cran.r-project.org/web/checks/check_results_vcfppR.html
+2. verbose message of vcfcomp

src/vcf-reader.cpp

@@ -121,22 +121,23 @@ using namespace std;
 //' }
 class vcfreader {
    public:
-    vcfreader(const std::string& vcffile) {
+    vcfreader(const std::string& vcffile) : fin(vcffile) {
         br.open(vcffile);
-        var.init(br.header);
+        var.initHeader(br.header);
     }
 
-    vcfreader(const std::string& vcffile, const std::string& region) {
+    vcfreader(const std::string& vcffile, const std::string& region) : fin(vcffile) {
         br.open(vcffile);
         if (!region.empty()) br.setRegion(region);
-        var.init(br.header);
+        var.initHeader(br.header);
     }
 
-    vcfreader(const std::string& vcffile, const std::string& region, const std::string& samples) {
+    vcfreader(const std::string& vcffile, const std::string& region, const std::string& samples)
+        : fin(vcffile) {
         br.open(vcffile);
         if (!samples.empty()) br.setSamples(samples);
         if (!region.empty()) br.setRegion(region);
-        var.init(br.header);
+        var.initHeader(br.header);
     }
 
     ~vcfreader() {}
@@ -259,11 +260,16 @@ class vcfreader {
     // WRITE
     inline void output(const std::string& vcffile) {
         bw.open(vcffile);
-        bw.initalHeader(br.header);
+        bw.copyHeader(fin);
+        var.resetHeader(bw.header);
         writable = true;
     }
-    inline void write() { bw.writeRecord(var); }
-    inline void close() { bw.close(); }
+    inline void write() {
+        if (writable) bw.writeRecord(var);
+    }
+    inline void close() {
+        if (writable) bw.close();
+    }
 
     inline void setCHR(std::string s) { var.setCHR(s.c_str()); }
     inline void setID(std::string s) { var.setID(s.c_str()); }
@@ -318,6 +324,7 @@ class vcfreader {
 
    private:
     bool writable = false;
+    const std::string fin;
     vcfpp::BcfReader br;
     vcfpp::BcfRecord var;
     vcfpp::BcfWriter bw;

src/vcfpp.h

@@ -2,7 +2,7 @@
  * @file        https://github.com/Zilong-Li/vcfpp/vcfpp.h
  * @author      Zilong Li
  * @email       [email protected]
- * @version     v0.3.6
+ * @version     v0.3.7
  * @breif       a single C++ file for manipulating VCF
  * Copyright (C) 2022-2023.The use of this code is governed by the LICENSE file.
  ******************************************************************************/
@@ -118,6 +118,15 @@ inline bool isEndWith(std::string const & s, std::string const & e)
     }
 }
 
+// determinate the mode for read/write the compressed/uncompressed VCF/BCF
+inline std::string getMode(std::string const & fname, std::string mode = "r")
+{
+    if(isEndWith(fname, "bcf.gz")) mode += "b";
+    if(isEndWith(fname, "bcf")) mode += "bu";
+    if(isEndWith(fname, "vcf.gz")) mode += "z";
+    return mode;
+}
+
 // string split by separator
 inline std::vector<std::string> split_string(const std::string & s, const std::string & separators)
 {
@@ -416,7 +425,7 @@ class BcfRecord
     /// constructor with a given BcfHeader object
     BcfRecord(BcfHeader & h)
     {
-        init(h);
+        initHeader(h);
     }
 
     ~BcfRecord()
@@ -425,8 +434,8 @@ class BcfRecord
         if(hdr_d) bcf_hdr_destroy(hdr_d);
     }
 
-    /// initilize a BcfRecord object by pointing to another BcfHeader object
-    void init(BcfHeader & h)
+    /// initilize the header associated with BcfRecord object by pointing to another BcfHeader object
+    void initHeader(BcfHeader & h)
     {
         header = &h;
         if(!header->hdr) throw std::runtime_error("please initial header first\n");
@@ -435,6 +444,12 @@ class BcfRecord
         gtPhase.resize(nsamples, 0);
     }
 
+    /// reset the header associated with BcfRecord object by pointing to another BcfHeader object
+    void resetHeader(BcfHeader & h)
+    {
+        header = &h;
+    }
+
     /** @brief stream out the variant */
     friend std::ostream & operator<<(std::ostream & out, const BcfRecord & v)
     {
@@ -1369,15 +1384,6 @@ class BcfReader
         SamplesName = header.getSamples();
     }
 
-    /// return if the file is opened successfully
-    bool isOpen() const
-    {
-        if(fp != NULL)
-            return true;
-        else
-            return false;
-    }
-
     /** @brief set the number of threads to use */
     inline int setThreads(int n)
     {
@@ -1562,10 +1568,7 @@ class BcfWriter
      */
     void open(const std::string & fname)
     {
-        std::string mode{"w"};
-        if(isEndWith(fname, "bcf.gz")) mode += "b";
-        if(isEndWith(fname, "bcf")) mode += "bu";
-        if(isEndWith(fname, "vcf.gz")) mode += "z";
+        auto mode = getMode(fname, "w");
         fp = std::shared_ptr<htsFile>(hts_open(fname.c_str(), mode.c_str()), htsFile_close());
     }
 
@@ -1604,6 +1607,15 @@ class BcfWriter
         hp = &h;
     }
 
+    /// copy header of given VCF
+    void copyHeader(const std::string & vcffile)
+    {
+        htsFile * fp2 = hts_open(vcffile.c_str(), "r");
+        header.hdr = bcf_hdr_read(fp2);
+        hts_close(fp2);
+        initalHeader(header);
+    }
+
     /// copy a string to a vcf line
     void writeLine(const std::string & vcfline)
     {

tests/testthat/test-modify-vcf.R

@@ -24,7 +24,6 @@ test_that("modify the genotypes", {
   g2 <- br$genotypes(F)
   expect_false(isTRUE(all.equal(g0, g2)))
   expect_identical(g1, g2)
-  br$close()
   ##  the original vcf can not be modified
   br <- vcfreader$new(outvcf)
   br$variant() ## get a variant record
@@ -54,7 +53,6 @@ test_that("modify item in FORMAT", {
   ds <- gp[2,] + gp[3,] * 2
   ## will prompt uses a message if `output` is not called 
   br$addFORMAT("DS", "1", "Float", "Diploid dosage")
-  br$addINFO("INFO", "1", "Float", "INFO score of imputation")
   ## now open another file for output 
   newvcf <- paste0(tempfile(), ".vcf.gz")
   br$output(newvcf)