zzh
/
minimap2
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

r369: updated example with the latest API

This commit is contained in:
Heng Li 2017-09-14 22:44:10 -04:00
parent 75ff7ceec5
commit 14b853499f
2 changed files with 36 additions and 38 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

33
example.c
View File

@ -11,7 +11,13 @@ KSEQ_INIT(gzFile, gzread)
int main(int argc, char *argv[]) int main(int argc, char *argv[])
{ {
mm_idxopt_t iopt;
mm_mapopt_t mopt;
int n_threads = 3;
mm_verbose = 2; // disable message output to stderr mm_verbose = 2; // disable message output to stderr
mm_set_opt(0, &iopt, &mopt);
mopt.flag |= MM_F_CIGAR; // perform alignment
if (argc < 3) { if (argc < 3) {
fprintf(stderr, "Usage: minimap2-lite <target.fa> <query.fa>\n"); fprintf(stderr, "Usage: minimap2-lite <target.fa> <query.fa>\n");
@ -23,24 +29,17 @@ int main(int argc, char *argv[])
assert(f); assert(f);
kseq_t *ks = kseq_init(f); kseq_t *ks = kseq_init(f);
// create index for target; we are creating one index for all target sequence // open index reader
int n_threads = 4, w = 10, k = 15, is_hpc = 0; mm_idx_reader_t *r = mm_idx_reader_open(argv[1], &iopt, 0);
mm_idx_t *mi = mm_idx_build(argv[1], w, k, is_hpc, n_threads); mm_idx_t *mi;
assert(mi); while ((mi = mm_idx_reader_read(r, n_threads)) != 0) { // traverse each part of the index
mm_mapopt_update(&mopt, mi); // this sets the maximum minimizer occurrence; TODO: set a better default in mm_mapopt_init()!
// mapping
mm_mapopt_t opt;
mm_mapopt_init(&opt); // initialize mapping parameters
mm_mapopt_update(&opt, mi); // this sets the maximum minimizer occurrence; TODO: set a better default in mm_mapopt_init()!
opt.flag |= MM_F_CIGAR; // perform alignment
mm_tbuf_t *tbuf = mm_tbuf_init(); // thread buffer; for multi-threading, allocate one tbuf for each thread mm_tbuf_t *tbuf = mm_tbuf_init(); // thread buffer; for multi-threading, allocate one tbuf for each thread
while (kseq_read(ks) >= 0) { // each kseq_read() call reads one query sequence while (kseq_read(ks) >= 0) { // each kseq_read() call reads one query sequence
mm_reg1_t *reg; mm_reg1_t *reg;
int j, i, n_reg; int j, i, n_reg;
// get all hits for the query reg = mm_map(mi, ks->seq.l, ks->seq.s, &n_reg, tbuf, &mopt, 0); // get all hits for the query
reg = mm_map(mi, ks->seq.l, ks->seq.s, &n_reg, tbuf, &opt, 0); for (j = 0; j < n_reg; ++j) { // traverse hits and print them out
// traverse hits and print them out
for (j = 0; j < n_reg; ++j) {
mm_reg1_t *r = &reg[j]; mm_reg1_t *r = &reg[j];
assert(r->p); // with MM_F_CIGAR, this should not be NULL assert(r->p); // with MM_F_CIGAR, this should not be NULL
printf("%s\t%d\t%d\t%d\t%c\t", ks->name.s, ks->seq.l, r->qs, r->qe, "+-"[r->rev]); printf("%s\t%d\t%d\t%d\t%c\t", ks->name.s, ks->seq.l, r->qs, r->qe, "+-"[r->rev]);
@ -54,10 +53,10 @@ int main(int argc, char *argv[])
free(reg); free(reg);
} }
mm_tbuf_destroy(tbuf); mm_tbuf_destroy(tbuf);
// deallocate index and close the query file
mm_idx_destroy(mi); mm_idx_destroy(mi);
kseq_destroy(ks); }
mm_idx_reader_close(r); // close the index reader
kseq_destroy(ks); // close the query file
gzclose(f); gzclose(f);
return 0; return 0;
} }

7
main.c
View File

@ -6,7 +6,7 @@
#include "mmpriv.h" #include "mmpriv.h"
#include "getopt.h" #include "getopt.h"
#define MM_VERSION "2.1.1-r368-dirty" #define MM_VERSION "2.1.1-r369-dirty"
#ifdef __linux__ #ifdef __linux__
#include <sys/resource.h> #include <sys/resource.h>
@ -67,6 +67,7 @@ int main(int argc, char *argv[])
char *fnw = 0, *rg = 0, *s; char *fnw = 0, *rg = 0, *s;
FILE *fp_help = stderr; FILE *fp_help = stderr;
mm_idx_reader_t *idx_rdr; mm_idx_reader_t *idx_rdr;
mm_idx_t *mi;
liftrlimit(); liftrlimit();
mm_realtime0 = realtime(); mm_realtime0 = realtime();
@ -210,9 +211,7 @@ int main(int argc, char *argv[])
} }
if (opt.flag & MM_F_OUT_SAM) if (opt.flag & MM_F_OUT_SAM)
mm_write_sam_hdr_no_SQ(rg, MM_VERSION, argc, argv); mm_write_sam_hdr_no_SQ(rg, MM_VERSION, argc, argv);
for (;;) { while ((mi = mm_idx_reader_read(idx_rdr, n_threads)) != 0) {
mm_idx_t *mi;
if ((mi = mm_idx_reader_read(idx_rdr, n_threads)) == 0) break;
if (mm_verbose >= 2 && idx_rdr->n_parts > 1 && (opt.flag&MM_F_OUT_SAM) && !(opt.flag&MM_F_NO_SAM_SQ)) if (mm_verbose >= 2 && idx_rdr->n_parts > 1 && (opt.flag&MM_F_OUT_SAM) && !(opt.flag&MM_F_NO_SAM_SQ))
fprintf(stderr, "[WARNING] \033[1;31mSAM output is malformated due to internal @SQ lines. Please add option --no-sam-sq or filter afterwards.\033[0m\n"); fprintf(stderr, "[WARNING] \033[1;31mSAM output is malformated due to internal @SQ lines. Please add option --no-sam-sq or filter afterwards.\033[0m\n");
if (mm_verbose >= 3) if (mm_verbose >= 3)