Parallel bgzip with cram support

Makefiles/Makefile.include

@@ -72,4 +72,4 @@ CFLAGS ?= $(OPTFLAG) -pipe -Wall
 COMPFLAGS = $(CFLAGS) $(USER_WARNINGS) -I$(INCLUDE_PATH) $(CURRENT_DIR_INCLUDE) $(USER_INCLUDES) $(USE_KNET) $(USE_ZLIB) -D_FILE_OFFSET_BITS=64 -D__STDC_LIMIT_MACROS $(USER_COMPILE_VARS)
 
 # default installation directory
-INSTALLDIR?=/usr/local/bin
+INSTALLDIR?=/usr/local/lib

Makefiles/Makefile.test

@@ -38,7 +38,7 @@ all debug: $(EXE)
 
 # dependencies for executables
 $(EXE) : $(LIBRARY) $(OBJECTS)
-	$(CXX) $(COMPFLAGS) -o  $@ $(OBJECTS) $(LIBRARY) $(LDFLAGS) -lm $(ZLIB_LIB) $(UNAME_LIBS)
+	$(CXX) $(COMPFLAGS) -o  $@ $(OBJECTS) $(LIBRARY) $(LDFLAGS) -lm -lhts $(ZLIB_LIB) $(UNAME_LIBS)
 
 $(OBJECTS): $(TOOLHDR) $(LIBHDR) | $(OBJDIR)
 

README.md

@@ -5,6 +5,10 @@ On debian type systems (including Ubuntu), add the following packages if they ar
 
     sudo apt-get install g++ libssl-dev zlib1g-dev
 
+The following libraries are dependencies that need to be installed on the build system.
+
+[htslib 1.9](https://github.com/samtools/htslib/tree/1.9)
+
 Building
 --------
 
@@ -23,7 +27,6 @@ Under the main statgen repository, there are:
 - `glf` - library code for operating on glf files.
 - `include` - after compiling, the library headers are linked here
 - `Makefiles` - directory containing Makefiles that are used in the library and can be used for developing programs using the library
-- `samtools` - library code used from samtools
 
 After Compiling: `libStatGen.a`, `libStatGen_debug.a`, and `libStatGen_profile.a` are created at the top level.
 

general/BgzfFileType.cpp

@@ -27,27 +27,56 @@ bool BgzfFileType::ourRequireEofBlock = true;
 
 BgzfFileType::BgzfFileType(const char * filename, const char * mode)
 {
+    BgzfFileType::numThreads = 1;
+    char threadSpec[8];
+    char bgzfMode[8];
+    const char* thread_start = strchr(mode, '@');
+    int mode_len = 0;
+    if (thread_start != NULL){
+        mode_len = thread_start-mode;
+    } else {
+        mode_len = strnlen(mode, 7);
+    }
+
+    if(mode_len < 8){
+        strncpy(bgzfMode, mode, mode_len);
+        bgzfMode[mode_len] = '\0';
+    } else {
+        strncpy(bgzfMode, mode, 7);
+        bgzfMode[7] = '\0';
+    }
+    if (thread_start != NULL && thread_start[1] != '\0'){
+        // Advance past the @, then parse the number of threads
+        thread_start++;
+        strncpy(threadSpec, thread_start, 7);
+        threadSpec[7] = '\0';
+        numThreads = strtol(threadSpec, NULL, 10);
+        // If we can't parse the number, revert to one thread
+        if (numThreads == 0){ numThreads = 1; }
+    }
     // If the file is for write and is '-', then write to stdout.
-    if(((mode[0] == 'w') || (mode[0] == 'W')) && 
+    if(((mode[0] == 'w') || (mode[0] == 'W')) &&
        (strcmp(filename, "-") == 0))
     {
         // Write to stdout.
-        bgzfHandle = bgzf_dopen(fileno(stdout), mode);
+        bgzfHandle = bgzf_dopen(fileno(stdout), bgzfMode);
     }
-    else if(((mode[0] == 'r') || (mode[0] == 'R')) && 
+    else if(((mode[0] == 'r') || (mode[0] == 'R')) &&
        (strcmp(filename, "-") == 0))
     {
         // read from stdin
-        bgzfHandle = bgzf_dopen(fileno(stdin), mode);
+        bgzfHandle = bgzf_dopen(fileno(stdin), bgzfMode);
     }
     else
     {
-        bgzfHandle = bgzf_open(filename, mode);
+        bgzfHandle = bgzf_open(filename, bgzfMode);
     }
 
     myStartPos = 0;
     if (bgzfHandle != NULL)
     {
+        //Only do multithreaded IO if more than one thread is used.
+        if(numThreads > 1){ bgzf_mt(bgzfHandle, numThreads, 256); }
         // Check to see if the file is being opened for read, if the eof block
         // is required, and if it is, if it is there.
         if ((mode[0] == 'r' || mode[0] == 'R') && (strcmp(filename, "-") != 0)

general/BgzfFileType.h

@@ -29,12 +29,14 @@ class BgzfFileType : public FileType
 public:
     BgzfFileType()
     {
+        numThreads = 1;
         bgzfHandle = NULL;
         myEOF = false;
     }
 
     virtual ~BgzfFileType()
     {
+        bgzf_close(bgzfHandle);
         bgzfHandle = NULL;
     }
 
@@ -110,7 +112,7 @@ class BgzfFileType : public FileType
         }
         else if((bytesRead != (int)size) & (bytesRead >= 0))
         {
-            // Less then the requested size was read 
+            // Less then the requested size was read
             // and an error was not returned (bgzf_read returns -1 on error).
             myEOF = true;
         }
@@ -173,6 +175,9 @@ class BgzfFileType : public FileType
     // at the end of the file.  If the block is required, but not on the file,
     // the constructor fails to open the file.
     static bool ourRequireEofBlock;
+
+    // The number of threads to use when using multithreaded BGZF
+    int64_t numThreads;
 };
 
 #endif

general/InputFile.cpp

@@ -32,9 +32,21 @@ InputFile::InputFile(const char * filename, const char * mode,
     myAttemptRecovery = false;
     myFileTypePtr = NULL;
     myBufferIndex = 0;
+    myWriteIndex = 0;
     myCurrentBufferSize = 0;
     myAllocatedBufferSize = DEFAULT_BUFFER_SIZE;
-    myFileBuffer = new char[myAllocatedBufferSize];
+    if(strchr(mode, 'r') || strchr(mode, 'R')){
+        myFileBuffer = new char[myAllocatedBufferSize];
+    } else {
+        myFileBuffer = NULL;
+    }
+    if(strchr(mode, 'w') || strchr(mode, 'W') ||
+       strchr(mode, 'a') || strchr(mode, 'A')){
+        myWriteBuffer = new char[myAllocatedBufferSize];
+        memset(myWriteBuffer, '\0', myAllocatedBufferSize);
+    } else {
+        myWriteBuffer = NULL;
+    }
     myFileName.clear();
 
     openFile(filename, mode, compressionMode);
@@ -55,13 +67,13 @@ int InputFile::readTilChar(const std::string& stopChars, std::string& stringRef)
         {
             return(-1);
         }
-        
+
         // Try to find the character in the stopChars.
         pos = stopChars.find(charRead);
 
         if(pos == std::string::npos)
         {
-            // Didn't find a stop character and it is not an EOF, 
+            // Didn't find a stop character and it is not an EOF,
             // so add it to the string.
             stringRef += charRead;
         }
@@ -84,7 +96,7 @@ int InputFile::readTilChar(const std::string& stopChars)
         {
             return(-1);
         }
-        
+
         // Try to find the character in the stopChars.
         pos = stopChars.find(charRead);
     }
@@ -163,7 +175,7 @@ bool InputFile::openFile(const char * filename, const char * mode,
 {
     //
     // if recovering, we don't want to issue big readaheads, since
-    // that interferes with the decompression - we only want to 
+    // that interferes with the decompression - we only want to
     // decompress one at a time, and handle the exceptions immediately
     // rather than at some indeterminate point in time.
     //
@@ -180,7 +192,7 @@ bool InputFile::openFile(const char * filename, const char * mode,
         // Check if reading from stdin.
         if((strcmp(filename, "-") == 0) || (strcmp(filename, "-.gz") == 0))
         {
-            // Reading from stdin, open it based on the 
+            // Reading from stdin, open it based on the
             // compression mode.
             openFileUsingMode(filename, mode, compressionMode);
         }
@@ -219,7 +231,7 @@ bool InputFile::openFile(const char * filename, const char * mode,
                 // Read the file to see if it a gzip file.
                 GzipHeader gzipHeader;
                 bool isGzip = gzipHeader.readHeader(file);
-                
+
                 // The file header has been read, so close the file, so it can
                 // be re-opened as the correct type.
                 file.close();
@@ -392,6 +404,11 @@ InputFile::~InputFile()
         delete[] myFileBuffer;
         myFileBuffer = NULL;
     }
+    if(myWriteBuffer != NULL)
+    {
+        delete[] myWriteBuffer;
+        myWriteBuffer = NULL;
+    }
 }