From 716713685ce6379e70ad9e2f0ab065065421b865 Mon Sep 17 00:00:00 2001 From: Heng Li Date: Wed, 25 Mar 2020 19:59:17 -0400 Subject: [PATCH] renamed files for code cleanup --- Assembly.cpp | 2 +- Hash_Table.h | 2 +- Makefile | 25 +++++++------- Trio.cpp | 2 +- hist.cpp | 2 +- yak-count.cpp => htab.cpp | 69 +++++++++++++++++++++++++++++++++++---- yak.h => htab.h | 13 ++------ main.cpp | 2 +- sketch.cpp | 2 +- yak-sys.cpp => sys.cpp | 2 +- yak-bbf.cpp | 42 ------------------------ 11 files changed, 83 insertions(+), 80 deletions(-) rename yak-count.cpp => htab.cpp (90%) rename yak.h => htab.h (90%) rename yak-sys.cpp => sys.cpp (97%) delete mode 100644 yak-bbf.cpp diff --git a/Assembly.cpp b/Assembly.cpp index 7c1b527..f0a3aa6 100644 --- a/Assembly.cpp +++ b/Assembly.cpp @@ -9,7 +9,7 @@ #include "POA.h" #include "Correct.h" #include "Output.h" -#include "yak.h" +#include "htab.h" Total_Count_Table TCB; Total_Pos_Table PCB; diff --git a/Hash_Table.h b/Hash_Table.h index c4ae1c8..7090fdb 100644 --- a/Hash_Table.h +++ b/Hash_Table.h @@ -2,7 +2,7 @@ #define __HASHTABLE__ #include "khashl.h" #include "kmer.h" -#include "yak.h" +#include "htab.h" KHASHL_MAP_INIT(static inline, Count_Table, ha_ct, uint64_t, int, kh_hash_dummy, kh_eq_generic) KHASHL_MAP_INIT(static inline, Pos_Table, ha_pt, uint64_t, uint64_t, kh_hash_dummy, kh_eq_generic) diff --git a/Makefile b/Makefile index d411cc9..6925f23 100644 --- a/Makefile +++ b/Makefile @@ -4,7 +4,7 @@ CPPFLAGS= INCLUDES= OBJS= Output.o CommandLines.o Process_Read.o Assembly.o Hash_Table.o \ POA.o Correct.o Levenshtein_distance.o Overlaps.o Trio.o kthread.o \ - yak-bbf.o yak-count.o hist.o sketch.o yak-sys.o + htab.o hist.o sketch.o sys.o EXE= hifiasm LIBS= -lz -lpthread -lm @@ -33,30 +33,29 @@ depend: # DO NOT DELETE Assembly.o: Assembly.h CommandLines.h Process_Read.h kseq.h Overlaps.h kvec.h -Assembly.o: kdq.h kmer.h Hash_Table.h khashl.h yak.h POA.h Correct.h +Assembly.o: kdq.h kmer.h Hash_Table.h khashl.h htab.h POA.h Correct.h Assembly.o: Levenshtein_distance.h Output.h CommandLines.o: CommandLines.h ketopt.h Correct.o: Correct.h Hash_Table.h khashl.h kmer.h Process_Read.h kseq.h -Correct.o: Overlaps.h kvec.h kdq.h CommandLines.h yak.h +Correct.o: Overlaps.h kvec.h kdq.h CommandLines.h htab.h Correct.o: Levenshtein_distance.h POA.h Assembly.h Hash_Table.o: Hash_Table.h khashl.h kmer.h Process_Read.h kseq.h Overlaps.h -Hash_Table.o: kvec.h kdq.h CommandLines.h yak.h Correct.h +Hash_Table.o: kvec.h kdq.h CommandLines.h htab.h Correct.h Hash_Table.o: Levenshtein_distance.h POA.h ksort.h Levenshtein_distance.o: Levenshtein_distance.h Output.o: Output.h CommandLines.h Overlaps.o: Overlaps.h kvec.h kdq.h ksort.h Process_Read.h kseq.h -Overlaps.o: CommandLines.h Hash_Table.h khashl.h kmer.h yak.h Correct.h +Overlaps.o: CommandLines.h Hash_Table.h khashl.h kmer.h htab.h Correct.h Overlaps.o: Levenshtein_distance.h POA.h POA.o: POA.h Hash_Table.h khashl.h kmer.h Process_Read.h kseq.h Overlaps.h -POA.o: kvec.h kdq.h CommandLines.h yak.h Correct.h Levenshtein_distance.h +POA.o: kvec.h kdq.h CommandLines.h htab.h Correct.h Levenshtein_distance.h Process_Read.o: Process_Read.h kseq.h Overlaps.h kvec.h kdq.h CommandLines.h Trio.o: khashl.h kthread.h Process_Read.h kseq.h Overlaps.h kvec.h kdq.h -Trio.o: CommandLines.h yak.h -hist.o: yak.h +Trio.o: CommandLines.h htab.h +hist.o: htab.h +htab.o: kthread.h khashl.h kseq.h htab.h CommandLines.h kthread.o: kthread.h main.o: CommandLines.h Process_Read.h kseq.h Overlaps.h kvec.h kdq.h -main.o: Assembly.h Levenshtein_distance.h yak.h -sketch.o: kvec.h yak.h -yak-bbf.o: yak.h -yak-count.o: CommandLines.h yak.h khashl.h kthread.h kseq.h -yak-sys.o: yak.h +main.o: Assembly.h Levenshtein_distance.h htab.h +sketch.o: kvec.h htab.h +sys.o: htab.h diff --git a/Trio.cpp b/Trio.cpp index 27b3dc8..17eb05e 100644 --- a/Trio.cpp +++ b/Trio.cpp @@ -6,7 +6,7 @@ #include "khashl.h" // hash table #include "kthread.h" #include "Process_Read.h" -#include "yak.h" +#include "htab.h" #include "CommandLines.h" #define YAK_MAX_KMER 31 diff --git a/hist.cpp b/hist.cpp index 728990d..5021c4b 100644 --- a/hist.cpp +++ b/hist.cpp @@ -1,5 +1,5 @@ #include -#include "yak.h" +#include "htab.h" static void yak_hist_line(int c, int x, int exceed, int64_t cnt) { diff --git a/yak-count.cpp b/htab.cpp similarity index 90% rename from yak-count.cpp rename to htab.cpp index 4eeb488..b7595e1 100644 --- a/yak-count.cpp +++ b/htab.cpp @@ -7,7 +7,7 @@ #include "kthread.h" #include "khashl.h" #include "kseq.h" -#include "yak.h" +#include "htab.h" #include "CommandLines.h" #define YAK_COUNTER_BITS 12 @@ -57,6 +57,54 @@ void yak_copt_init(yak_copt_t *o) o->chunk_size = 10000000; } +/************************ + * Blocked bloom filter * + ************************/ + +typedef struct { + int n_shift, n_hashes; + uint8_t *b; +} yak_bf_t; + +yak_bf_t *yak_bf_init(int n_shift, int n_hashes) +{ + yak_bf_t *b; + void *ptr = 0; + if (n_shift + YAK_BLK_SHIFT > 64 || n_shift < YAK_BLK_SHIFT) return 0; + CALLOC(b, 1); + b->n_shift = n_shift; + b->n_hashes = n_hashes; + posix_memalign(&ptr, 1<<(YAK_BLK_SHIFT-3), 1ULL<<(n_shift-3)); + b->b = (uint8_t*)ptr; + bzero(b->b, 1ULL<<(n_shift-3)); + return b; +} + +void yak_bf_destroy(yak_bf_t *b) +{ + if (b == 0) return; + free(b->b); free(b); +} + +int yak_bf_insert(yak_bf_t *b, uint64_t hash) +{ + int x = b->n_shift - YAK_BLK_SHIFT; + uint64_t y = hash & ((1ULL<> x & YAK_BLK_MASK; + int h2 = hash >> b->n_shift & YAK_BLK_MASK; + uint8_t *p = &b->b[y<<(YAK_BLK_SHIFT-3)]; + int i, z = h1, cnt = 0; + if ((h2&31) == 0) h2 = (h2 + 1) & YAK_BLK_MASK; // otherwise we may repeatedly use a few bits + for (i = 0; i < b->n_hashes; z = (z + h2) & YAK_BLK_MASK) { + uint8_t *q = &p[z>>3], u; + u = 1<<(z&7); + cnt += !!(*q & u); + *q |= u; + ++i; + } + return cnt; +} + /******************** * Count hash table * ********************/ @@ -306,7 +354,7 @@ typedef struct { ha_mz1_t *b; } ch_buf_t; -static inline void ch_insert_buf(ch_buf_t *buf, int p, uint64_t y) // insert a k-mer $y to a linear buffer +static inline void ct_insert_buf(ch_buf_t *buf, int p, uint64_t y) // insert a k-mer $y to a linear buffer { int pre = y & ((1<> 1 | (uint64_t)(1 - (c&1)) << shift; x[3] = x[3] >> 1 | (uint64_t)(1 - (c>>1)) << shift; if (++l >= k) - ch_insert_buf(buf, p, yak_hash_long(x)); + ct_insert_buf(buf, p, yak_hash_long(x)); } else l = 0, x[0] = x[1] = x[2] = x[3] = 0; // if there is an "N", restart } } @@ -358,7 +406,7 @@ static void count_seq_buf_HPC(ch_buf_t *buf, int k, int p, int len, const char * x[2] = x[2] >> 1 | (uint64_t)(1 - (c&1)) << shift; x[3] = x[3] >> 1 | (uint64_t)(1 - (c>>1)) << shift; if (++l >= k) - ch_insert_buf(buf, p, yak_hash_long(x)); + ct_insert_buf(buf, p, yak_hash_long(x)); last = c; } } else l = 0, last = -1, x[0] = x[1] = x[2] = x[3] = 0; // if there is an "N", restart @@ -375,6 +423,8 @@ typedef struct { // global data structure for kt_pipeline() const yak_copt_t *opt; const void *flt_tab; int create_new, is_store; + uint64_t batch_offset; + uint64_t n_base; kseq_t *ks; ha_ct_t *ct; ha_pt_t *pt; @@ -405,7 +455,7 @@ static void worker_for_mz(void *data, long i, int tid) st_data_t *s = (st_data_t*)data; ha_mz1_v *b = &s->mz_buf[tid]; s->mz_buf[tid].n = 0; - ha_sketch(s->seq[i], s->len[i], s->p->opt->w, s->p->opt->k, 0, s->p->opt->is_HPC, b, s->p->flt_tab); + ha_sketch(s->seq[i], s->len[i], s->p->opt->w, s->p->opt->k, s->p->batch_offset + i, s->p->opt->is_HPC, b, s->p->flt_tab); s->mz[i].n = s->mz[i].m = b->n; MALLOC(s->mz[i].a, b->n); memcpy(s->mz[i].a, b->a, b->n * sizeof(ha_mz1_t)); @@ -421,7 +471,11 @@ static void *worker_count(void *data, int step, void *in) // callback for kt_pip s->p = p; while ((ret = kseq_read(p->ks)) >= 0) { int l = p->ks->seq.l; - if (l < p->opt->k) continue; + if (p->batch_offset + s->n_seq >= (1<<28) - 1) { + fprintf(stderr, "ERROR: this implementation supports no more than %d reads\n", (1<<28) - 1); + exit(1); + } + p->n_base += l; if (s->n_seq == s->m_seq) { s->m_seq = s->m_seq < 16? 16 : s->m_seq + (s->m_seq>>1); REALLOC(s->len, s->m_seq); @@ -473,7 +527,7 @@ static void *worker_count(void *data, int step, void *in) // callback for kt_pip } else { for (i = 0; i < s->n_seq; ++i) for (j = 0; j < s->mz[i].n; ++j) - ch_insert_buf(s->buf, p->opt->pre, s->mz[i].a[j].x); + ct_insert_buf(s->buf, p->opt->pre, s->mz[i].a[j].x); } for (i = 0; i < s->n_seq; ++i) { free(s->mz[i].a); @@ -494,6 +548,7 @@ static void *worker_count(void *data, int step, void *in) // callback for kt_pip free(s->buf[i].a); } p->ct->tot += n_ins; + p->batch_offset += s->n_seq; free(s->buf); fprintf(stderr, "[M::%s::%.3f*%.2f] processed %d sequences; %ld distinct k-mers in the hash table\n", __func__, yak_realtime(), yak_cputime() / yak_realtime(), s->n_seq, (long)p->ct->tot); diff --git a/yak.h b/htab.h similarity index 90% rename from yak.h rename to htab.h index 52373af..dbb844b 100644 --- a/yak.h +++ b/htab.h @@ -1,5 +1,5 @@ -#ifndef __YAK_H__ -#define __YAK_H__ +#ifndef __HA_HTAB_H__ +#define __HA_HTAB_H__ #define __STDC_LIMIT_MACROS #include @@ -17,11 +17,6 @@ typedef struct { typedef struct { uint32_t n, m; ha_mz1_t *a; } ha_mz1_v; -typedef struct { - int n_shift, n_hashes; - uint8_t *b; -} yak_bf_t; - extern const unsigned char seq_nt4_table[256]; int ha_ft_isflt(const void *hh, uint64_t y); @@ -37,10 +32,6 @@ long yak_peakrss(void); void ha_sketch(const char *str, int len, int w, int k, uint32_t rid, int is_hpc, ha_mz1_v *p, const void *hf); int yak_analyze_count(int n_cnt, const int64_t *cnt, int *peak_het); -yak_bf_t *yak_bf_init(int n_shift, int n_hashes); -void yak_bf_destroy(yak_bf_t *b); -int yak_bf_insert(yak_bf_t *b, uint64_t hash); - static inline uint64_t yak_hash64(uint64_t key, uint64_t mask) // invertible integer hash function { key = (~key + (key << 21)) & mask; // key = (key << 21) - key - 1; diff --git a/main.cpp b/main.cpp index d5f8995..1e71192 100644 --- a/main.cpp +++ b/main.cpp @@ -4,7 +4,7 @@ #include "Process_Read.h" #include "Assembly.h" #include "Levenshtein_distance.h" -#include "yak.h" +#include "htab.h" int main(int argc, char *argv[]) { diff --git a/sketch.cpp b/sketch.cpp index 58c64f4..c3546ed 100644 --- a/sketch.cpp +++ b/sketch.cpp @@ -3,7 +3,7 @@ #include #include #include "kvec.h" -#include "yak.h" +#include "htab.h" typedef struct { // a simplified version of kdq int front, count; diff --git a/yak-sys.cpp b/sys.cpp similarity index 97% rename from yak-sys.cpp rename to sys.cpp index d01c68c..fed8a2a 100644 --- a/yak-sys.cpp +++ b/sys.cpp @@ -1,6 +1,6 @@ #include #include -#include "yak.h" +#include "htab.h" int yak_verbose = 3; diff --git a/yak-bbf.cpp b/yak-bbf.cpp deleted file mode 100644 index 742b80a..0000000 --- a/yak-bbf.cpp +++ /dev/null @@ -1,42 +0,0 @@ -#include -#include -#include "yak.h" - -yak_bf_t *yak_bf_init(int n_shift, int n_hashes) -{ - yak_bf_t *b; - void *ptr = 0; - if (n_shift + YAK_BLK_SHIFT > 64 || n_shift < YAK_BLK_SHIFT) return 0; - CALLOC(b, 1); - b->n_shift = n_shift; - b->n_hashes = n_hashes; - posix_memalign(&ptr, 1<<(YAK_BLK_SHIFT-3), 1ULL<<(n_shift-3)); - b->b = (uint8_t*)ptr; - bzero(b->b, 1ULL<<(n_shift-3)); - return b; -} - -void yak_bf_destroy(yak_bf_t *b) -{ - if (b == 0) return; - free(b->b); free(b); -} - -int yak_bf_insert(yak_bf_t *b, uint64_t hash) -{ - int x = b->n_shift - YAK_BLK_SHIFT; - uint64_t y = hash & ((1ULL<> x & YAK_BLK_MASK; - int h2 = hash >> b->n_shift & YAK_BLK_MASK; - uint8_t *p = &b->b[y<<(YAK_BLK_SHIFT-3)]; - int i, z = h1, cnt = 0; - if ((h2&31) == 0) h2 = (h2 + 1) & YAK_BLK_MASK; // otherwise we may repeatedly use a few bits - for (i = 0; i < b->n_hashes; z = (z + h2) & YAK_BLK_MASK) { - uint8_t *q = &p[z>>3], u; - u = 1<<(z&7); - cnt += !!(*q & u); - *q |= u; - ++i; - } - return cnt; -}