from wifionICE

Author: pep04706

include/bloom.h

@@ -54,6 +54,21 @@ typedef struct _bfilter {
   unsigned char *filter; /**< filter sequence*/
 } Bfilter;
 
+
+/**
+ * @brief global parameters needed for the processed reads
+ *
+ * */
+typedef struct _bloom_par {
+  int kmersize;  /**< kmer size (number of elements)*/
+  int hashNum;  /**< number of hash functions used to construct the filter*/
+  int kmersizeBytes;  /**< Bytes needed to store the kmer (4bases ~ 1byte) */
+  int halfsizeBytes;  /**< half size in bytes(needed to decide whether
+                           to store a kmer or its reverse complement) */
+  int hangingBases;  /**< number of hanging bases that don't complete a byte*/
+  int hasOverhead;  /**< kmer has overhead when kmersize % 4!=0 */
+} Bloom_par;
+
 /**
  * @brief stores a processed kmer (2 bits pro nucleotide)
  *
@@ -79,6 +94,8 @@ void free_Bfilter(Bfilter *ptr_bf);
 
 Procs_kmer *init_procs(int kmersize, int hashNum);
 
+Procs_kmer *init_procs(int kmersize, int hashNum);
+
 void free_procs(Procs_kmer *procs);
 
 double score_read_in_filter(unsigned char *read, int L, Procs_kmer *procs,

include/defines.h

@@ -29,52 +29,52 @@
 #ifndef DEFINES_H_
 #define DEFINES_H_
 
-// General 
+// General
 #define B_LEN 131072  /**< buffer size */
 #define MAX_FILENAME 300 /**< Maximum # chars in a filename */
-#define bool short  /**< define a bool type */
+#define bool int16_t  /**< define a bool type */
 #define true 1   /**< assign true to 1 */
 #define false 0  /**< assign false to 0 */
 
 #ifndef max
-   #define max( a, b ) ( ((a) > (b)) ? (a) : (b) ) /**< max function */
+  #define max(a, b) (((a) > (b)) ? (a) : (b)) /**< max function */
 #endif
 
 #ifndef min
-   #define min( a, b ) ( ((a) < (b)) ? (a) : (b) ) /**< min function */
+  #define min(a, b) (((a) < (b)) ? (a) : (b)) /**< min function */
 #endif
 
 #ifndef mem_usageMB
-   #define mem_usageMB()  fprintf(stderr, \
+  #define mem_usageMB()  fprintf(stderr, \
          "- Current allocated memory: %ld MB.\n", \
          alloc_mem >> 20)  /**< returns allocated memory in MB */
 #endif
 
 #ifndef mem_usage
-   #define mem_usage()  fprintf(stderr, \
+  #define mem_usage()  fprintf(stderr, \
          "- Current allocated memory: %ld Bytes.\n", \
          alloc_mem) /**< returns allocated memory in Bytes */
 #endif
 
 
 // Q_report, S_report
-#define DEFAULT_MINQ 27   /**< Minimum quality threshold */ 
-#define DEFAULT_NTILES 96 /**< Default number of tiles */ 
-#define DEFAULT_NQ 46 /**< Default number of different quality values */ 
-#define ZEROQ 33 /**< ASCII code of lowest quality value (!) */ 
-#define N_ACGT 5 /**< Number of different nucleotides in the fq file */ 
+#define DEFAULT_MINQ 27   /**< Minimum quality threshold */
+#define DEFAULT_NTILES 96 /**< Default number of tiles */
+#define DEFAULT_NQ 46 /**< Default number of different quality values */
+#define ZEROQ 33 /**< ASCII code of lowest quality value (!) */
+#define N_ACGT 5 /**< Number of different nucleotides in the fq file */
 #define MAX_RCOMMAND  4000 /**< Maximum # chars in R command*/
 
 
 // Fasta files
 #define FA_ENTRY_BUF 20 /**< buffer for fasta entries*/
 
-// Tree 
+// Tree
 #define T_ACGT 4 /**< Number of children per node in tree*/
 #define NPOOL_1D 1048576 /**< Number of Node structs allocated in inner dim */
 #define NPOOL_2D 16  /**< Number of *Node allocated in outer dim */
 #define MAX_FASZ_TREE 1e7 /**< Maximum fasta size for constructing a tree.
-                               DECIDE A SENSIBLE SIZE! */
+                               DECIDE A SENSIBLE SIZE */
 // BloomFilter
 #define BITSPERCHAR 8  /**< number of bits in a char */
 #define BASESPERCHAR 4  /**< number of nucleotides that can fit in a char */
@@ -90,17 +90,17 @@
 #define STRIP 3     /**< Looks for the largest N-free sequence */
 // trimQ only
 #define FRAC  3     /**< Discards a read if it contains > percent lowQ bases*/
-#define ENDSFRAC 4  /**< trims at the ends and discards a read if the 
+#define ENDSFRAC 4  /**< trims at the ends and discards a read if the
                       remaining part has more than > percent lowQ bases */
-#define GLOBAL 5    /**< Trims a fixed # bases from e left and right*/    
+#define GLOBAL 5    /**< Trims a fixed # bases from e left and right*/
 
 #define TREE 1   /**< Use a tree to look for contaminations*/
 #define SA 2     /**< Use a suffix array to look for contaminations*/
 #define BLOOM 3  /**< Use a bloom filter to look for contaminations*/
 
 #define ERROR 1000 /**< Encodes an error when reading in trimN, trimQ, method
                      options in trimFilter */
-#define DEFAULT_MINL 25 /**< Default minimum length under which we discard 
+#define DEFAULT_MINL 25 /**< Default minimum length under which we discard
                           the reads */
 
 // Classification of filters
@@ -110,8 +110,8 @@
 #define NNNN 3 /**<  N's presence filter */
 #define GOOD 4 /**<  Good reads */
 
-//Number of filters
-#define NFILTERS 4 /**< total number of filters */ 
+// Number of filters
+#define NFILTERS 4 /**< total number of filters */
 
 #endif // endif DEFINES_H_
 

include/init_trimFilter.h

@@ -30,9 +30,9 @@
 #ifndef INIT_TRIMFILTER_H_
 #define INIT_TRIMFILTER_H_
 
+#include <uintstd.h>
 #include "defines.h"
 
-
 /**
  * @ brief adapter struct
  * @ note UNFINISHED!

src/bloom.c

@@ -139,6 +139,30 @@ void free_Bfilter(Bfilter * ptr_bf) {
   alloc_mem -= ptr_bf -> bfsizeBytes;
 }
 
+/** 
+ * @brief initializes bloom parameters
+ * @param kmersize number of elements of the kmer
+ * @param hashNum number of hash functions to be computed
+ * @return Bloom_par structure
+ *
+ * */
+Bloom_par init_bloom_par(int kmersize, int hashNum) {
+  Bloom_par res; 
+  res.kmersizeBytes = kmersize / 4;
+  res.halfsizeBytes = kmersize / 8;
+  res.hangingBases = 0;
+  res.hasOverhead = 0;
+  res.hashNum = hashNum;
+  if (kmersize % 8 != 0) {
+      res.halfsizeBytes++;
+      if ((res.hangingBases = kmersize % 4) > 0) {
+            res.kmersizeBytes++;
+            res.hasOverhead = 1;
+      }
+  }
+  return res; 
+}
+
 /**
  * @brief initializes a Procs structure, given the kmersize and the
  *        number of hash functions
@@ -448,7 +472,7 @@ bool insert_and_fetch(Bfilter *ptr_bf, Procs_kmer* procs) {
   int i = 0;
   uint64_t modValue;
   // iterates through hashed values adding it to the filter
-  for (i = 0; i < procs -> hashNum; i++) {
+  for (i = 0; i < ptr_bf -> hashNum; i++) {
      modValue = (procs -> hashValues[i]) % (ptr_bf -> bfsizeBits);
      result &= ((__sync_fetch_and_or(&(ptr_bf->filter[modValue/BITSPERCHAR]),
            bitMask[modValue % BITSPERCHAR]))>>(modValue % BITSPERCHAR)) & 1;