diff --git a/Assembly.h b/Assembly.h index 44c7ef9..ef7caba 100644 --- a/Assembly.h +++ b/Assembly.h @@ -9,6 +9,8 @@ #define Get_Cigar_Length(RECORD) (RECORD>>2) void *ha_gen_flt_tab(const hifiasm_opt_t *asm_opt); +void *ha_gen_mzidx(const hifiasm_opt_t *asm_opt, const void *flt_tab); + void Counting_multiple_thr(); void Build_hash_table_multiple_thr(); void Overlap_calculate_multipe_thr(); diff --git a/CommandLines.cpp b/CommandLines.cpp index b402e30..3a6213a 100644 --- a/CommandLines.cpp +++ b/CommandLines.cpp @@ -65,6 +65,7 @@ void init_opt(hifiasm_opt_t* asm_opt) asm_opt->mat_index = NULL; asm_opt->thread_num = 1; asm_opt->k_mer_length = 41; + asm_opt->mz_win = 41; asm_opt->bf_shift = 37; asm_opt->high_factor = 10.0f; asm_opt->no_HPC = 0; diff --git a/CommandLines.h b/CommandLines.h index 17426dd..2dfbbd5 100644 --- a/CommandLines.h +++ b/CommandLines.h @@ -15,6 +15,7 @@ typedef struct { char* mat_index; int thread_num; int k_mer_length; + int mz_win; int bf_shift; float high_factor; int no_HPC; diff --git a/main.cpp b/main.cpp index 4ffd8a7..cb5a2b8 100644 --- a/main.cpp +++ b/main.cpp @@ -17,6 +17,7 @@ int main(int argc, char *argv[]) yak_reset_realtime(); flt_tab = ha_gen_flt_tab(&asm_opt); + ha_gen_mzidx(&asm_opt, flt_tab); ha_hf_destroy(flt_tab); if (0) { Correct_Reads(asm_opt.number_of_round); diff --git a/yak-count.cpp b/yak-count.cpp index 7c2fc47..c7062b3 100644 --- a/yak-count.cpp +++ b/yak-count.cpp @@ -496,3 +496,24 @@ void ha_hf_destroy(void *h) { yak_hh_destroy((yak_hh_t*)h); } + +void *ha_gen_mzidx(const hifiasm_opt_t *asm_opt, const void *flt_tab) +{ + int64_t cnt[YAK_N_COUNTS]; + int peak_hom, peak_het; + yak_copt_t opt; + yak_ch_t *h; + yak_copt_init(&opt); + opt.is_HPC = !asm_opt->no_HPC; + opt.k = asm_opt->k_mer_length; + opt.w = asm_opt->mz_win; + opt.n_thread = asm_opt->thread_num; + opt.bf_shift = asm_opt->bf_shift; + h = yak_count_file(&opt, asm_opt->num_reads, asm_opt->read_file_names, flt_tab); + yak_ch_hist(h, cnt, opt.n_thread); + fprintf(stderr, "[M::%s] count[%d] = %ld (for sanity check)\n", __func__, YAK_MAX_COUNT, (long)cnt[YAK_MAX_COUNT]); + peak_hom = yak_analyze_count(YAK_N_COUNTS, cnt, &peak_het); + if (peak_hom > 0) fprintf(stderr, "[M::%s] peak_hom: %d; peak_het: %d\n", __func__, peak_hom, peak_het); + yak_ch_destroy(h); + return 0; +}