20 const int32_t nodeIdx,
const int32_t depth,
const int32_t offset,
21 const int32_t y0,
float64_t*
const*
const W,
const int32_t K )
25 float64_t*
const W_kiy = & W[ depth ][ offset + y0 ];
26 POIMTrie*
const node = &TreeMem[ nodeIdx ];
29 if( depth < degree-1 )
31 const int32_t offset1 = offset * NUM_SYMS;
32 for( sym = 0; sym < NUM_SYMS; ++sym )
35 const int32_t childIdx = node->children[ sym ];
36 if( childIdx != NO_CHILD )
38 W_kiy[ sym ] = TreeMem[ childIdx ].weight;
42 const int32_t y1 = ( y0 + sym ) * NUM_SYMS;
43 POIMs_extract_W_helper( childIdx, depth+1, offset1, y1, W, K );
51 for( sym = 0; sym < NUM_SYMS; ++sym )
54 W_kiy[ sym ] = node->child_weights[ sym ];
61 float64_t*
const*
const W,
const int32_t K)
65 const int32_t N = length;
67 for( i = 0; i < N; ++i ) {
69 POIMs_extract_W_helper( trees[i], 0, i*NUM_SYMS, 0*NUM_SYMS, W, K );
75 const float64_t*
const distrib,
const int32_t i,
const int32_t nodeIdx,
76 int32_t left_tries_idx[4],
const int32_t depth, int32_t
const lastSym,
82 const float64_t*
const distribLeft = & distrib[ (i-1) * NUM_SYMS ];
83 POIMTrie*
const node = &TreeMem[ nodeIdx ];
94 const float64_t*
const distribRight = & distrib[ (i+depth) * NUM_SYMS ];
97 for( symRight = 0; symRight < NUM_SYMS; ++symRight ) {
98 const float64_t w1 = node->child_weights[ symRight ];
99 const float64_t pRight = distribRight[ symRight ];
107 for( symLeft = 0; symLeft < NUM_SYMS; ++symLeft )
109 if (left_tries_idx[symLeft] != NO_CHILD)
111 POIMTrie* nodeLeft = &TreeMem[left_tries_idx[symLeft]];
114 const float64_t w2 = nodeLeft->child_weights[ lastSym ];
115 const float64_t pLeft = distribLeft[ symLeft ];
135 const float64_t*
const distrib,
const int32_t i,
const int32_t nodeIdx,
139 ASSERT(0<=depth && depth<=degree-2)
142 const float64_t*
const distribLeft = & distrib[ (i-1) * NUM_SYMS ];
143 POIMTrie*
const node = &TreeMem[ nodeIdx ];
156 for( symRight = 0; symRight < NUM_SYMS; ++symRight )
158 const int32_t childIdx = node->
children[ symRight ];
159 if( childIdx != NO_CHILD )
161 if( depth < degree-2 )
163 int32_t new_left_tries_idx[4];
165 for( symLeft = 0; symLeft < NUM_SYMS; ++symLeft )
167 new_left_tries_idx[symLeft]=NO_CHILD;
169 if (left_tries_idx[symLeft] != NO_CHILD)
171 POIMTrie* nodeLeft = &TreeMem[left_tries_idx[symLeft]];
173 new_left_tries_idx[ symLeft ]=nodeLeft->children[ symRight ];
176 POIMs_calc_SLR_helper2( distrib, i, childIdx, new_left_tries_idx, depth+1, &auxS, &dummy, &auxR );
179 POIMs_calc_SLR_helper1( distrib, i, childIdx, left_tries_idx, depth+1, symRight, &auxS, &dummy, &auxR );
183 const float64_t*
const distribRight = & distrib[ (i+depth) * NUM_SYMS ];
184 const float64_t pRight = distribRight[ symRight ];
185 node->S += pRight * auxS;
186 node->R += pRight * auxR;
192 for( symLeft = 0; symLeft < NUM_SYMS; ++symLeft )
194 if (left_tries_idx[symLeft] != NO_CHILD)
196 const POIMTrie* nodeLeft = &TreeMem[left_tries_idx[symLeft]];
198 const float64_t pLeft = distribLeft[ symLeft ];
200 node->S += pLeft * nodeLeft->S;
201 node->L += pLeft * nodeLeft->L;
205 const float64_t*
const distribRight = & distrib[ (i+depth) * NUM_SYMS ];
208 if( depth < degree-2 )
210 for( symRight = 0; symRight < NUM_SYMS; ++symRight )
212 const int32_t childIdxLeft = nodeLeft->children[ symRight ];
213 if( childIdxLeft != NO_CHILD )
215 const POIMTrie*
const childLeft = &TreeMem[ childIdxLeft ];
216 auxS += distribRight[symRight] * childLeft->S;
222 for( symRight = 0; symRight < NUM_SYMS; ++symRight ) {
223 auxS += distribRight[symRight] * nodeLeft->child_weights[ symRight ];
226 node->S -= pLeft* auxS;
252 const int32_t N = length;
255 int32_t leftSubtrees[4];
256 for( symLeft = 0; symLeft < NUM_SYMS; ++symLeft )
257 leftSubtrees[ symLeft ] = NO_CHILD;
259 for(int32_t i = 0; i < N; ++i )
261 POIMs_calc_SLR_helper2( distrib, i, trees[i], leftSubtrees, 0, &dummy, &dummy, &dummy );
263 const POIMTrie*
const node = &TreeMem[ trees[i] ];
264 ASSERT(trees[i]!=NO_CHILD)
266 for(symLeft = 0; symLeft < NUM_SYMS; ++symLeft )
267 leftSubtrees[ symLeft ] = node->children[ symLeft ];
273 const int32_t parentIdx,
const int32_t sym,
const int32_t depth,
276 ASSERT(parentIdx!=NO_CHILD)
277 const POIMTrie*
const parent = &TreeMem[ parentIdx ];
278 if( depth < degree ) {
279 const int32_t nodeIdx = parent->children[ sym ];
280 const POIMTrie*
const node = &TreeMem[ nodeIdx ];
286 const float64_t w = parent->child_weights[ sym ];
295 float64_t*
const*
const poims,
const int32_t K,
const int32_t k,
296 const int32_t i,
const int32_t y,
const float64_t valW,
301 const int32_t nk = nofsKmers[ k ];
308 if( debug==0 || debug==2 )
310 poims[ k-1 ][ i*nk + y ] += valS - valW;
314 if( debug==0 || debug==3 )
317 for( r = 1; k+r <= K; ++r )
319 const int32_t nr = nofsKmers[ r ];
320 const int32_t nz = nofsKmers[ k+r ];
321 float64_t*
const poim = & poims[ k+r-1 ][ i*nz ];
323 for( j = 0; j < nr; ++j )
325 if( !( 0 <= z && z < nz ) ) {
326 SG_PRINT(
"k=%d, nk=%d, r=%d, nr=%d, nz=%d \n", k, nk, r, nr, nz )
327 SG_PRINT(
" j=%d, y=%d, z=%d \n", j, y, z )
330 poim[ z ] += valL - valW;
336 if( debug==0 || debug==4 )
339 for( l = 1; k+l <= K && l <= i; ++l )
341 const int32_t nl = nofsKmers[ l ];
342 const int32_t nz = nofsKmers[ k+l ];
343 float64_t*
const poim = & poims[ k+l-1 ][ (i-l)*nz ];
345 for( j = 0; j < nl; ++j ) {
347 poim[ z ] += valR - valW;
356 const int32_t nodeIdx,
const int32_t depth,
const int32_t i,
357 const int32_t y0,
float64_t*
const*
const poims,
const int32_t K,
362 POIMTrie*
const node = &TreeMem[ nodeIdx ];
364 if( depth < degree-1 )
368 for( sym = 0; sym < NUM_SYMS; ++sym )
370 const int32_t childIdx = node->
children[ sym ];
371 if( childIdx != NO_CHILD )
373 POIMTrie*
const child = &TreeMem[ childIdx ];
374 const int32_t y = y0 + sym;
375 const int32_t y1 = y * NUM_SYMS;
377 POIMs_add_SLR_helper2( poims, K, depth+1, i, y, w, child->S, child->L, child->R, debug );
378 POIMs_add_SLR_helper1( childIdx, depth+1, i, y1, poims, K, debug );
385 for( sym = 0; sym < NUM_SYMS; ++sym )
387 const int32_t childIdx = node->children[ sym ];
388 if( childIdx != NO_CHILD ) {
389 POIMTrie*
const child = &TreeMem[ childIdx ];
390 const int32_t y = y0 + sym;
392 POIMs_add_SLR_helper2( poims, K, depth+1, i, y, w, child->S, child->L, child->R, debug );
400 for( sym = 0; sym < NUM_SYMS; ++sym )
402 const float64_t w = node->child_weights[ sym ];
403 const int32_t y = y0 + sym;
404 POIMs_add_SLR_helper2( poims, K, depth+1, i, y, w, w, w, w, debug );
413 float64_t*
const*
const poims,
const int32_t K,
const int32_t debug)
418 const int32_t m = ( degree > K ) ? degree : K;
419 nofsKmers = SG_MALLOC(int32_t, m+1 );
423 for( k = 0; k < m+1; ++k ) {
428 const int32_t N = length;
430 for( i = 0; i < N; ++i ) {
431 POIMs_add_SLR_helper1( trees[i], 0, i, 0*NUM_SYMS, poims, K, debug );
Template class Trie implements a suffix trie, i.e. a tree in which all suffixes up to a certain lengt...
void POIMs_extract_W_helper(const int32_t nodeIdx, const int32_t depth, const int32_t offset, const int32_t y0, float64_t *const *const W, const int32_t K)
all of classes and functions are contained in the shogun namespace