|
0
|
1 #include <stdlib.h>
|
|
|
2 #include <stdio.h>
|
|
|
3 #include <assert.h>
|
|
|
4 #include "bam.h"
|
|
|
5 #include "ksort.h"
|
|
|
6
|
|
|
7 #define TV_GAP 2
|
|
|
8
|
|
|
9 typedef struct __freenode_t {
|
|
|
10 uint32_t level:28, cnt:4;
|
|
|
11 struct __freenode_t *next;
|
|
|
12 } freenode_t, *freenode_p;
|
|
|
13
|
|
|
14 #define freenode_lt(a,b) ((a)->cnt < (b)->cnt || ((a)->cnt == (b)->cnt && (a)->level < (b)->level))
|
|
|
15 KSORT_INIT(node, freenode_p, freenode_lt)
|
|
|
16
|
|
|
17 /* Memory pool, similar to the one in bam_pileup.c */
|
|
|
18 typedef struct {
|
|
|
19 int cnt, n, max;
|
|
|
20 freenode_t **buf;
|
|
|
21 } mempool_t;
|
|
|
22
|
|
|
23 static mempool_t *mp_init()
|
|
|
24 {
|
|
|
25 return (mempool_t*)calloc(1, sizeof(mempool_t));
|
|
|
26 }
|
|
|
27 static void mp_destroy(mempool_t *mp)
|
|
|
28 {
|
|
|
29 int k;
|
|
|
30 for (k = 0; k < mp->n; ++k) free(mp->buf[k]);
|
|
|
31 free(mp->buf); free(mp);
|
|
|
32 }
|
|
|
33 static inline freenode_t *mp_alloc(mempool_t *mp)
|
|
|
34 {
|
|
|
35 ++mp->cnt;
|
|
|
36 if (mp->n == 0) return (freenode_t*)calloc(1, sizeof(freenode_t));
|
|
|
37 else return mp->buf[--mp->n];
|
|
|
38 }
|
|
|
39 static inline void mp_free(mempool_t *mp, freenode_t *p)
|
|
|
40 {
|
|
|
41 --mp->cnt; p->next = 0; p->cnt = TV_GAP;
|
|
|
42 if (mp->n == mp->max) {
|
|
|
43 mp->max = mp->max? mp->max<<1 : 256;
|
|
|
44 mp->buf = (freenode_t**)realloc(mp->buf, sizeof(freenode_t*) * mp->max);
|
|
|
45 }
|
|
|
46 mp->buf[mp->n++] = p;
|
|
|
47 }
|
|
|
48
|
|
|
49 /* core part */
|
|
|
50 struct __bam_lplbuf_t {
|
|
|
51 int max, n_cur, n_pre;
|
|
|
52 int max_level, *cur_level, *pre_level;
|
|
|
53 mempool_t *mp;
|
|
|
54 freenode_t **aux, *head, *tail;
|
|
|
55 int n_nodes, m_aux;
|
|
|
56 bam_pileup_f func;
|
|
|
57 void *user_data;
|
|
|
58 bam_plbuf_t *plbuf;
|
|
|
59 };
|
|
|
60
|
|
|
61 void bam_lplbuf_reset(bam_lplbuf_t *buf)
|
|
|
62 {
|
|
|
63 freenode_t *p, *q;
|
|
|
64 bam_plbuf_reset(buf->plbuf);
|
|
|
65 for (p = buf->head; p->next;) {
|
|
|
66 q = p->next;
|
|
|
67 mp_free(buf->mp, p);
|
|
|
68 p = q;
|
|
|
69 }
|
|
|
70 buf->head = buf->tail;
|
|
|
71 buf->max_level = 0;
|
|
|
72 buf->n_cur = buf->n_pre = 0;
|
|
|
73 buf->n_nodes = 0;
|
|
|
74 }
|
|
|
75
|
|
|
76 static int tview_func(uint32_t tid, uint32_t pos, int n, const bam_pileup1_t *pl, void *data)
|
|
|
77 {
|
|
|
78 bam_lplbuf_t *tv = (bam_lplbuf_t*)data;
|
|
|
79 freenode_t *p;
|
|
|
80 int i, l, max_level;
|
|
|
81 // allocate memory if necessary
|
|
|
82 if (tv->max < n) { // enlarge
|
|
|
83 tv->max = n;
|
|
|
84 kroundup32(tv->max);
|
|
|
85 tv->cur_level = (int*)realloc(tv->cur_level, sizeof(int) * tv->max);
|
|
|
86 tv->pre_level = (int*)realloc(tv->pre_level, sizeof(int) * tv->max);
|
|
|
87 }
|
|
|
88 tv->n_cur = n;
|
|
|
89 // update cnt
|
|
|
90 for (p = tv->head; p->next; p = p->next)
|
|
|
91 if (p->cnt > 0) --p->cnt;
|
|
|
92 // calculate cur_level[]
|
|
|
93 max_level = 0;
|
|
|
94 for (i = l = 0; i < n; ++i) {
|
|
|
95 const bam_pileup1_t *p = pl + i;
|
|
|
96 if (p->is_head) {
|
|
|
97 if (tv->head->next && tv->head->cnt == 0) { // then take a free slot
|
|
|
98 freenode_t *p = tv->head->next;
|
|
|
99 tv->cur_level[i] = tv->head->level;
|
|
|
100 mp_free(tv->mp, tv->head);
|
|
|
101 tv->head = p;
|
|
|
102 --tv->n_nodes;
|
|
|
103 } else tv->cur_level[i] = ++tv->max_level;
|
|
|
104 } else {
|
|
|
105 tv->cur_level[i] = tv->pre_level[l++];
|
|
|
106 if (p->is_tail) { // then return a free slot
|
|
|
107 tv->tail->level = tv->cur_level[i];
|
|
|
108 tv->tail->next = mp_alloc(tv->mp);
|
|
|
109 tv->tail = tv->tail->next;
|
|
|
110 ++tv->n_nodes;
|
|
|
111 }
|
|
|
112 }
|
|
|
113 if (tv->cur_level[i] > max_level) max_level = tv->cur_level[i];
|
|
|
114 ((bam_pileup1_t*)p)->level = tv->cur_level[i];
|
|
|
115 }
|
|
|
116 assert(l == tv->n_pre);
|
|
|
117 tv->func(tid, pos, n, pl, tv->user_data);
|
|
|
118 // sort the linked list
|
|
|
119 if (tv->n_nodes) {
|
|
|
120 freenode_t *q;
|
|
|
121 if (tv->n_nodes + 1 > tv->m_aux) { // enlarge
|
|
|
122 tv->m_aux = tv->n_nodes + 1;
|
|
|
123 kroundup32(tv->m_aux);
|
|
|
124 tv->aux = (freenode_t**)realloc(tv->aux, sizeof(void*) * tv->m_aux);
|
|
|
125 }
|
|
|
126 for (p = tv->head, i = l = 0; p->next;) {
|
|
|
127 if (p->level > max_level) { // then discard this entry
|
|
|
128 q = p->next;
|
|
|
129 mp_free(tv->mp, p);
|
|
|
130 p = q;
|
|
|
131 } else {
|
|
|
132 tv->aux[i++] = p;
|
|
|
133 p = p->next;
|
|
|
134 }
|
|
|
135 }
|
|
|
136 tv->aux[i] = tv->tail; // add a proper tail for the loop below
|
|
|
137 tv->n_nodes = i;
|
|
|
138 if (tv->n_nodes) {
|
|
|
139 ks_introsort(node, tv->n_nodes, tv->aux);
|
|
|
140 for (i = 0; i < tv->n_nodes; ++i) tv->aux[i]->next = tv->aux[i+1];
|
|
|
141 tv->head = tv->aux[0];
|
|
|
142 } else tv->head = tv->tail;
|
|
|
143 }
|
|
|
144 // clean up
|
|
|
145 tv->max_level = max_level;
|
|
|
146 memcpy(tv->pre_level, tv->cur_level, tv->n_cur * 4);
|
|
|
147 // squeeze out terminated levels
|
|
|
148 for (i = l = 0; i < n; ++i) {
|
|
|
149 const bam_pileup1_t *p = pl + i;
|
|
|
150 if (!p->is_tail)
|
|
|
151 tv->pre_level[l++] = tv->pre_level[i];
|
|
|
152 }
|
|
|
153 tv->n_pre = l;
|
|
|
154 /*
|
|
|
155 fprintf(stderr, "%d\t", pos+1);
|
|
|
156 for (i = 0; i < n; ++i) {
|
|
|
157 const bam_pileup1_t *p = pl + i;
|
|
|
158 if (p->is_head) fprintf(stderr, "^");
|
|
|
159 if (p->is_tail) fprintf(stderr, "$");
|
|
|
160 fprintf(stderr, "%d,", p->level);
|
|
|
161 }
|
|
|
162 fprintf(stderr, "\n");
|
|
|
163 */
|
|
|
164 return 0;
|
|
|
165 }
|
|
|
166
|
|
|
167 bam_lplbuf_t *bam_lplbuf_init(bam_pileup_f func, void *data)
|
|
|
168 {
|
|
|
169 bam_lplbuf_t *tv;
|
|
|
170 tv = (bam_lplbuf_t*)calloc(1, sizeof(bam_lplbuf_t));
|
|
|
171 tv->mp = mp_init();
|
|
|
172 tv->head = tv->tail = mp_alloc(tv->mp);
|
|
|
173 tv->func = func;
|
|
|
174 tv->user_data = data;
|
|
|
175 tv->plbuf = bam_plbuf_init(tview_func, tv);
|
|
|
176 return (bam_lplbuf_t*)tv;
|
|
|
177 }
|
|
|
178
|
|
|
179 void bam_lplbuf_destroy(bam_lplbuf_t *tv)
|
|
|
180 {
|
|
|
181 freenode_t *p, *q;
|
|
|
182 free(tv->cur_level); free(tv->pre_level);
|
|
|
183 bam_plbuf_destroy(tv->plbuf);
|
|
|
184 free(tv->aux);
|
|
|
185 for (p = tv->head; p->next;) {
|
|
|
186 q = p->next;
|
|
|
187 mp_free(tv->mp, p); p = q;
|
|
|
188 }
|
|
|
189 mp_free(tv->mp, p);
|
|
|
190 assert(tv->mp->cnt == 0);
|
|
|
191 mp_destroy(tv->mp);
|
|
|
192 free(tv);
|
|
|
193 }
|
|
|
194
|
|
|
195 int bam_lplbuf_push(const bam1_t *b, bam_lplbuf_t *tv)
|
|
|
196 {
|
|
|
197 return bam_plbuf_push(b, tv->plbuf);
|
|
|
198 }
|
