diff --git a/.vscode/ipch/4eee59091dc293a3/PROCESS_READ.ipch b/.vscode/ipch/4eee59091dc293a3/PROCESS_READ.ipch index 1b643f6..8c78979 100644 Binary files a/.vscode/ipch/4eee59091dc293a3/PROCESS_READ.ipch and b/.vscode/ipch/4eee59091dc293a3/PROCESS_READ.ipch differ diff --git a/.vscode/ipch/756c62c8d5e4771d/COMMANDLINES.ipch b/.vscode/ipch/756c62c8d5e4771d/COMMANDLINES.ipch index 25872b9..943f0e4 100644 Binary files a/.vscode/ipch/756c62c8d5e4771d/COMMANDLINES.ipch and b/.vscode/ipch/756c62c8d5e4771d/COMMANDLINES.ipch differ diff --git a/.vscode/ipch/b49c1e8b711f08aa/HASH_TABLE.ipch b/.vscode/ipch/b49c1e8b711f08aa/HASH_TABLE.ipch index 1b88f27..b48acd6 100644 Binary files a/.vscode/ipch/b49c1e8b711f08aa/HASH_TABLE.ipch and b/.vscode/ipch/b49c1e8b711f08aa/HASH_TABLE.ipch differ diff --git a/.vscode/ipch/c0e71cfe49f0fe81/ASSEMBLY.ipch b/.vscode/ipch/c0e71cfe49f0fe81/ASSEMBLY.ipch index a9c0d8c..1b0ec93 100644 Binary files a/.vscode/ipch/c0e71cfe49f0fe81/ASSEMBLY.ipch and b/.vscode/ipch/c0e71cfe49f0fe81/ASSEMBLY.ipch differ diff --git a/.vscode/ipch/e64bc513afbf8204/mmap_address.bin b/.vscode/ipch/e64bc513afbf8204/mmap_address.bin new file mode 100644 index 0000000..862b842 Binary files /dev/null and b/.vscode/ipch/e64bc513afbf8204/mmap_address.bin differ diff --git a/.vscode/ipch/f2b98741f94c3e10/MAIN.ipch b/.vscode/ipch/f2b98741f94c3e10/MAIN.ipch index 238f6ea..ea20ef0 100644 Binary files a/.vscode/ipch/f2b98741f94c3e10/MAIN.ipch and b/.vscode/ipch/f2b98741f94c3e10/MAIN.ipch differ diff --git a/Assembly.cpp b/Assembly.cpp index f0f441e..abc8b48 100644 --- a/Assembly.cpp +++ b/Assembly.cpp @@ -396,6 +396,19 @@ void Build_hash_table_multiple_thr() fprintf(stdout, "%-30s%18.2f\n\n", "Build hash table time:", Get_T() - start_time); ///destory_Total_Count_Table(&TCB); + if (write_index_to_disk) + { + write_Total_Pos_Table(&PCB, read_file_name); + destory_Total_Pos_Table(&PCB); + load_Total_Pos_Table(&PCB, read_file_name); + + write_All_reads(&R_INF, read_file_name); + destory_All_reads(&R_INF); + load_All_reads(&R_INF, read_file_name); + } + + + } diff --git a/CommandLines.cpp b/CommandLines.cpp index c6d19df..dc70ddc 100644 --- a/CommandLines.cpp +++ b/CommandLines.cpp @@ -11,6 +11,9 @@ int thread_num = 1; int k_mer_length = 40; int k_mer_min_freq = 9; int k_mer_max_freq = 66; +int load_index_from_disk = 0; +int write_index_to_disk = 0; + double Get_T(void) { @@ -41,12 +44,14 @@ int CommandLine_process (int argc, char *argv[]) ketopt_t opt = KETOPT_INIT; int i, c; - while ((c = ketopt(&opt, argc, argv, 1, "ht:o:q:k:", longopts)) >= 0) { + while ((c = ketopt(&opt, argc, argv, 1, "ht:o:q:k:lw", longopts)) >= 0) { if (c == 100 || c == 'h') Print_H(); else if (c == 103 || c == 't') thread_num = atoi(opt.arg); else if (c == 102 || c == 'o') output_file_name = opt.arg; else if (c == 101 || c == 'q') read_file_name = opt.arg; else if (c == 'k') k_mer_length = atoi(opt.arg); + else if (c == 'l') load_index_from_disk = 1; + else if (c == 'w') write_index_to_disk = 1; else if (c == '?') printf("unknown opt: -%c\n", opt.opt? opt.opt : ':'); else if (c == ':') printf("missing arg: -%c\n", opt.opt? opt.opt : ':'); } diff --git a/CommandLines.h b/CommandLines.h index 2c3fe58..0f85670 100644 --- a/CommandLines.h +++ b/CommandLines.h @@ -10,7 +10,8 @@ extern int thread_num; extern int k_mer_length; extern int k_mer_min_freq; extern int k_mer_max_freq; - +extern int load_index_from_disk; +extern int write_index_to_disk; int CommandLine_process (int argc, char *argv[]); diff --git a/Hash_Table.cpp b/Hash_Table.cpp index 6e7ed70..48490ba 100644 --- a/Hash_Table.cpp +++ b/Hash_Table.cpp @@ -88,6 +88,105 @@ void destory_Total_Count_Table(Total_Count_Table* TCB) } +void destory_Total_Pos_Table(Total_Pos_Table* TCB) +{ + free(TCB->k_mer_index); + free(TCB->sub_h_lock); + free(TCB->pos); + + int i; + for (i = 0; i < TCB->size; i++) + { + kh_destroy(POS64, TCB->sub_h[i]); + } + free(TCB->sub_h); +} + + +void write_Total_Pos_Table(Total_Pos_Table* TCB, char* read_file_name) +{ + fprintf(stdout, "Writing index to disk ...... \n"); + char* index_name = (char*)malloc(strlen(read_file_name)+5); + sprintf(index_name, "%s.idx", read_file_name); + FILE* fp = fopen(index_name, "w"); + fwrite(&TCB->prefix_bits, sizeof(TCB->prefix_bits), 1, fp); + fwrite(&TCB->suffix_bits, sizeof(TCB->suffix_bits), 1, fp); + fwrite(&TCB->suffix_mode, sizeof(TCB->suffix_mode), 1, fp); + fwrite(&TCB->size, sizeof(TCB->size), 1, fp); + fwrite(&TCB->useful_k_mer, sizeof(TCB->useful_k_mer), 1, fp); + fwrite(&TCB->total_occ, sizeof(TCB->total_occ), 1, fp); + fwrite(TCB->k_mer_index, sizeof(uint64_t), TCB->useful_k_mer+1, fp); + fwrite(TCB->pos, sizeof(k_mer_pos), TCB->total_occ, fp); + + + int i; + for (i = 0; i < TCB->size; i++) + { + kh_write(POS64, TCB->sub_h[i], fp); + } + + free(index_name); + fclose(fp); + fprintf(stdout, "Index has been written.\n"); +} + + +void load_Total_Pos_Table(Total_Pos_Table* TCB, char* read_file_name) +{ + fprintf(stdout, "Loading index to disk ...... \n"); + char* index_name = (char*)malloc(strlen(read_file_name)+5); + sprintf(index_name, "%s.idx", read_file_name); + FILE* fp = fopen(index_name, "r"); + + fread(&TCB->prefix_bits, sizeof(TCB->prefix_bits), 1, fp); + fread(&TCB->suffix_bits, sizeof(TCB->suffix_bits), 1, fp); + fread(&TCB->suffix_mode, sizeof(TCB->suffix_mode), 1, fp); + fread(&TCB->size, sizeof(TCB->size), 1, fp); + fread(&TCB->useful_k_mer, sizeof(TCB->useful_k_mer), 1, fp); + fread(&TCB->total_occ, sizeof(TCB->total_occ), 1, fp); + + + if (TCB->useful_k_mer+1) + { + TCB->k_mer_index = (uint64_t*)malloc(sizeof(uint64_t)*(TCB->useful_k_mer+1)); + fread(TCB->k_mer_index, sizeof(uint64_t), TCB->useful_k_mer+1, fp); + } + else + { + TCB->k_mer_index = NULL; + } + + + if (TCB->total_occ) + { + TCB->pos = (k_mer_pos*)malloc(sizeof(k_mer_pos)*TCB->total_occ); + fread(TCB->pos, sizeof(k_mer_pos), TCB->total_occ, fp); + } + else + { + TCB->pos = NULL; + } + + TCB->sub_h_lock = (Hash_table_spin_lock*)malloc(sizeof(Hash_table_spin_lock)*TCB->size); + memset(TCB->sub_h_lock, 0, sizeof(Hash_table_spin_lock)*TCB->size); + + TCB->sub_h = (Pos_Table**)malloc(sizeof(Pos_Table*)*TCB->size); + + int i; + for (i = 0; i < TCB->size; i++) + { + init_Pos_Table(&(TCB->sub_h[i])); + TCB->sub_h_lock[i].lock = 0; + kh_load(POS64, TCB->sub_h[i], fp); + } + + free(index_name); + fclose(fp); + fprintf(stdout, "Index has been loaded.\n"); + +} + + void Traverse_Counting_Table(Total_Count_Table* TCB, Total_Pos_Table* PCB, int k_mer_min_freq, int k_mer_max_freq) { int i; @@ -98,7 +197,7 @@ void Traverse_Counting_Table(Total_Count_Table* TCB, Total_Pos_Table* PCB, int k PCB->useful_k_mer = 0; PCB->total_occ = 0; - init_Total_Pos_Table(PCB, TCB); + ///init_Total_Pos_Table(PCB, TCB); khint_t t; ///这就是个迭代器 int absent; diff --git a/Hash_Table.h b/Hash_Table.h index a3074f7..30b5b6c 100644 --- a/Hash_Table.h +++ b/Hash_Table.h @@ -292,6 +292,11 @@ void destory_Total_Count_Table(Total_Count_Table* TCB); void init_Count_Table(Count_Table** table); void init_Pos_Table(Count_Table** pre_table, Pos_Table** table); +void destory_Total_Pos_Table(Total_Pos_Table* TCB); +void write_Total_Pos_Table(Total_Pos_Table* TCB, char* read_file_name); +void load_Total_Pos_Table(Total_Pos_Table* TCB, char* read_file_name); + + void Traverse_Counting_Table(Total_Count_Table* TCB, Total_Pos_Table* PCB, int k_mer_min_freq, int k_mer_max_freq); diff --git a/Process_Read.cpp b/Process_Read.cpp index e8ae813..f70ae6d 100644 --- a/Process_Read.cpp +++ b/Process_Read.cpp @@ -41,6 +41,130 @@ void init_All_reads(All_reads* r) } +void destory_All_reads(All_reads* r) +{ + uint64_t i = 0; + for (i = 0; i < r->total_reads; i++) + { + if (r->N_site[i] != NULL) + { + free(r->N_site[i]); + } + } + free(r->N_site); + free(r->read); + free(r->name); + free(r->name_index); +} + + +void write_All_reads(All_reads* r, char* read_file_name) +{ + fprintf(stdout, "Writing reads to disk ...... \n"); + char* index_name = (char*)malloc(strlen(read_file_name)+5); + sprintf(index_name, "%s.bin", read_file_name); + FILE* fp = fopen(index_name, "w"); + fwrite(&r->index_size, sizeof(r->index_size), 1, fp); + fwrite(&r->name_index_size, sizeof(r->name_index_size), 1, fp); + fwrite(&r->total_reads, sizeof(r->total_reads), 1, fp); + fwrite(&r->total_reads_bases, sizeof(r->total_reads_bases), 1, fp); + fwrite(&r->total_name_length, sizeof(r->total_name_length), 1, fp); + + uint64_t i = 0; + uint64_t zero = 0; + for (i = 0; i < r->total_reads; i++) + { + if (r->N_site[i] != NULL) + { + ///这个实际上是N的个数 + fwrite(&r->N_site[i][0], sizeof(r->N_site[i][0]), 1, fp); + if (r->N_site[i][0]) + { + ///r->N_site[i]这实际是个长为r->N_site[i][0]+1 + ///这里从r->N_site[i] + 1写入了r->N_site[i][0]个元素 + fwrite(r->N_site[i]+1, sizeof(r->N_site[i][0]), r->N_site[i][0], fp); + } + } + else + { + fwrite(&zero, sizeof(zero), 1, fp); + } + + + + } + + fwrite(r->read, sizeof(uint8_t), (r->total_reads_bases/4 + r->total_reads + 5), fp); + fwrite(r->name, sizeof(char), r->total_name_length, fp); + fwrite(r->index, sizeof(uint64_t), r->index_size, fp); + fwrite(r->name_index, sizeof(uint64_t), r->name_index_size, fp); + + + free(index_name); + fclose(fp); + fprintf(stdout, "Reads has been written.\n"); +} + + + +void load_All_reads(All_reads* r, char* read_file_name) +{ + fprintf(stdout, "Loading reads to disk ...... \n"); + char* index_name = (char*)malloc(strlen(read_file_name)+5); + sprintf(index_name, "%s.bin", read_file_name); + FILE* fp = fopen(index_name, "r"); + fread(&r->index_size, sizeof(r->index_size), 1, fp); + fread(&r->name_index_size, sizeof(r->name_index_size), 1, fp); + fread(&r->total_reads, sizeof(r->total_reads), 1, fp); + fread(&r->total_reads_bases, sizeof(r->total_reads_bases), 1, fp); + fread(&r->total_name_length, sizeof(r->total_name_length), 1, fp); + + uint64_t i = 0; + uint64_t zero = 0; + r->N_site = (uint64_t**)malloc(sizeof(uint64_t*)*r->total_reads); + for (i = 0; i < r->total_reads; i++) + { + + fread(&zero, sizeof(zero), 1, fp); + + if (zero) + { + + r->N_site[i] = (uint64_t*)malloc(sizeof(uint64_t)*(zero + 1)); + r->N_site[i][0] = zero; + if (r->N_site[i][0]) + { + ///r->N_site[i]这实际是个长为r->N_site[i][0]+1 + ///这里从r->N_site[i] + 1写入了r->N_site[i][0]个元素 + fread(r->N_site[i]+1, sizeof(r->N_site[i][0]), r->N_site[i][0], fp); + } + } + else + { + r->N_site[i] = NULL; + } + + } + + r->read = (uint8_t*)malloc(sizeof(uint8_t)*(r->total_reads_bases/4 + r->total_reads + 5)); + fread(r->read, sizeof(uint8_t), (r->total_reads_bases/4 + r->total_reads + 5), fp); + + r->name = (char*)malloc(sizeof(char)*r->total_name_length); + fread(r->name, sizeof(char), r->total_name_length, fp); + + r->index = (uint64_t*)malloc(sizeof(uint64_t)*r->index_size); + fread(r->index, sizeof(uint64_t), r->index_size, fp); + + r->name_index = (uint64_t*)malloc(sizeof(uint64_t)*r->name_index_size); + fread(r->name_index, sizeof(uint64_t), r->name_index_size, fp); + + + free(index_name); + fclose(fp); + fprintf(stdout, "Reads has been loaded.\n"); +} + + inline void insert_read(All_reads* r, kstring_t* read, kstring_t* name) { @@ -68,7 +192,8 @@ void malloc_All_reads(All_reads* r) ///必须加r->total_reads r->read = (uint8_t*)malloc(sizeof(uint8_t)*(r->total_reads_bases/4 + r->total_reads + 5)); r->name = (char*)malloc(sizeof(char)*r->total_name_length); - r->N_site = (uint64_t**)malloc(sizeof(uint64_t*)*r->total_reads); + r->N_site = (uint64_t**)calloc(r->total_reads, sizeof(uint64_t*)); + } void init_UC_Read(UC_Read* r) diff --git a/Process_Read.h b/Process_Read.h index 9a537bb..22ea009 100644 --- a/Process_Read.h +++ b/Process_Read.h @@ -110,6 +110,9 @@ void compress_base(uint8_t* dest, char* src, uint64_t src_l, uint64_t** N_site_l void init_UC_Read(UC_Read* r); void recover_UC_Read(UC_Read* r, All_reads* R_INF, uint64_t ID); void reverse_complement(char* pattern, uint64_t length); +void write_All_reads(All_reads* r, char* read_file_name); +void load_All_reads(All_reads* r, char* read_file_name); +void destory_All_reads(All_reads* r); void Counting_block(); void destory_R_buffer_block(R_buffer_block* curr_sub_block); diff --git a/khash.h b/khash.h index f75f347..25b3883 100644 --- a/khash.h +++ b/khash.h @@ -128,6 +128,8 @@ int main() { #include #include #include +#include + /* compiler specific configuration */ @@ -206,7 +208,9 @@ static const double __ac_HASH_UPPER = 0.77; extern khint_t kh_get_##name(const kh_##name##_t *h, khkey_t key); \ extern int kh_resize_##name(kh_##name##_t *h, khint_t new_n_buckets); \ extern khint_t kh_put_##name(kh_##name##_t *h, khkey_t key, int *ret); \ - extern void kh_del_##name(kh_##name##_t *h, khint_t x); + extern void kh_del_##name(kh_##name##_t *h, khint_t x);\ + extern void kh_write_##name(kh_##name##_t *h, FILE* fp)\ + extern void kh_load_##name(kh_##name##_t *h, FILE* fp) #define __KHASH_IMPL(name, SCOPE, khkey_t, khval_t, kh_is_map, __hash_func, __hash_equal) \ SCOPE kh_##name##_t *kh_init_##name(void) { \ @@ -352,6 +356,35 @@ static const double __ac_HASH_UPPER = 0.77; __ac_set_isdel_true(h->flags, x); \ --h->size; \ } \ + } \ + SCOPE void kh_write_##name(kh_##name##_t *h, FILE* fp)\ + {\ + fwrite(&(h->n_buckets), sizeof(khint_t), 1, fp);\ + fwrite(&(h->size), sizeof(khint_t), 1, fp);\ + fwrite(&(h->n_occupied), sizeof(khint_t), 1, fp);\ + fwrite(&(h->upper_bound), sizeof(khint_t), 1, fp);\ + if (h->n_buckets)\ + {\ + fwrite(h->flags, sizeof(khint32_t), __ac_fsize(h->n_buckets), fp);\ + fwrite(h->keys, sizeof(khkey_t), h->n_buckets, fp);\ + fwrite(h->vals, sizeof(khval_t), h->n_buckets, fp);\ + }\ + } \ + SCOPE void kh_load_##name(kh_##name##_t *h, FILE* fp)\ + {\ + fread(&(h->n_buckets), sizeof(khint_t), 1, fp);\ + fread(&(h->size), sizeof(khint_t), 1, fp);\ + fread(&(h->n_occupied), sizeof(khint_t), 1, fp);\ + fread(&(h->upper_bound), sizeof(khint_t), 1, fp);\ + if (h->n_buckets)\ + {\ + h->flags = (khint32_t*)kmalloc(__ac_fsize(h->n_buckets) * sizeof(khint32_t));\ + fread(h->flags, sizeof(khint32_t), __ac_fsize(h->n_buckets), fp);\ + h->keys = (khkey_t*)kmalloc(sizeof(khkey_t)*h->n_buckets);\ + fread(h->keys, sizeof(khkey_t), h->n_buckets, fp);\ + h->vals = (khval_t*)kmalloc(sizeof(khval_t)*h->n_buckets);\ + fread(h->vals, sizeof(khval_t), h->n_buckets, fp);\ + }\ } #define KHASH_DECLARE(name, khkey_t, khval_t) \ @@ -624,4 +657,12 @@ typedef const char *kh_cstr_t; #define KHASH_MAP_INIT_STR(name, khval_t) \ KHASH_INIT(name, kh_cstr_t, khval_t, 1, kh_str_hash_func, kh_str_hash_equal) + + + +#define kh_write(name, h, fp) kh_write_##name(h, fp) + +#define kh_load(name, h, fp) kh_load_##name(h, fp) + + #endif /* __AC_KHASH_H */