« back to all changes in this revision
Viewing changes to kalign2_alignment_types.c
- browse files at revision 3
- compare with another revision
- download diff
- download tarball
- view history from revision 3
- Committer: Bazaar Package Importer
- Author(s): Andreas Tille, Charles Plessy, David Paleino, Andreas Tille
- Date: 2008-03-18 09:54:17 UTC
- mfrom: (2.1.1 gutsy)
- Revision ID: james.westby@ubuntu.com-20080318095417-cepxmgf2xsoxk84l
Tags: 2.03-2
[ Charles Plessy ]
* debian/control:
- Add Subversion repository URL to debian/control.
- Moved the Homepage: field out from the package's description.
- Enhances: t-coffee.
[ David Paleino ]
* debian/kalign.1 added - manpages are now statically generated
* debian/control:
- B-D updated (see above)
- added myself to Uploaders
- moved XS-Vcs-* fields to Vcs-*
* debian/rules:
- reflecting static build of manpages
- minor changes
* Updated to Standards-Version 3.7.3 (no changes needed)
[ Andreas Tille ]
* Added myself to Uploaders
* debian/control:
- Add Subversion repository URL to debian/control.
- Moved the Homepage: field out from the package's description.
- Enhances: t-coffee.
[ David Paleino ]
* debian/kalign.1 added - manpages are now statically generated
* debian/control:
- B-D updated (see above)
- added myself to Uploaders
- moved XS-Vcs-* fields to Vcs-*
* debian/rules:
- reflecting static build of manpages
- minor changes
* Updated to Standards-Version 3.7.3 (no changes needed)
[ Andreas Tille ]
* Added myself to Uploaders
- files added:
- Makefile.in
- files removed:
- Makefile
- files modified:
- README
added
removed
kalign2_alignment_types.c
1
/*
2
kalign2_alignment_types.c
3
4
Released under GPL - see the 'COPYING' file
5
6
Copyright (C) 2006 Timo Lassmann <timolassmann@gmail.com>
7
8
This program is free software; you can redistribute it and/or modify
9
it under the terms of the GNU General Public License as published by
10
the Free Software Foundation; either version 2 of the License, or
11
any later version.
12
13
This program is distributed in the hope that it will be useful,
14
but WITHOUT ANY WARRANTY; without even the implied warranty of
15
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16
GNU General Public License for more details.
17
18
You should have received a copy of the GNU General Public License
19
along with this program; if not, write to the Free Software
20
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21
22
Please send bug reports, comments etc. to:
23
timolassmann@gmail.com
24
*/
25
26
#include "kalign2.h"
27
28
int** default_alignment(struct alignment* aln,int* tree,float**submatrix, int** map)
29
{
30
struct dp_matrix *dp = 0;
31
int i,j,g,a,b,c;
32
int len_a;
33
int len_b;
34
float** profile = 0;
35
float* profa = 0;
36
float* profb = 0;
37
38
profile = malloc(sizeof(float*)*numprofiles);
39
for ( i = 0;i< numprofiles;i++){
40
profile[i] = 0;
41
}
42
43
map = malloc(sizeof(int*)*numprofiles);
44
for ( i = 0;i < numprofiles;i++){
45
map[i] = 0;
46
}
47
48
49
dp = dp_matrix_alloc(dp,511,511);
50
51
fprintf(stderr,"\nAlignment:\n");
52
53
//c = numseq;
54
for (i = 0; i < (numseq-1);i++){
55
a = tree[i*3];
56
b = tree[i*3+1];
57
c = tree[i*3+2];
58
fprintf(stderr,"\r%8.0f percent done",(float)(i) /(float)numseq * 100);
59
//fprintf(stderr,"Aligning:%d %d->%d %d %d\n",a,b,c,numseq,i);
60
len_a = aln->sl[a];
61
len_b = aln->sl[b];
62
dp = dp_matrix_realloc(dp,len_a,len_b);
63
64
map[c] = malloc(sizeof(int) * (len_a+len_b+2));
65
for (j = len_a+len_b+2;j--;){
66
map[c][j] = 0;
67
}
68
if (a < numseq){
69
profile[a] = make_profile(profile[a],aln->s[a],len_a,submatrix);
70
}
71
if (b < numseq){
72
profile[b] = make_profile(profile[b],aln->s[b],len_b,submatrix);
73
}
74
profa = profile[a]+64;
75
profb = profile[b]+64;
76
77
set_gap_penalties(profile[a],len_a,aln->nsip[b]);
78
set_gap_penalties(profile[b],len_b,aln->nsip[a]);
79
if(aln->nsip[a] == 1){
80
if(aln->nsip[b] == 1){
81
map[c] = ss_dyn(submatrix,map[c],dp,aln->s[a],aln->s[b],len_a,len_b);
82
}else{
83
map[c] = ps_dyn(map[c],dp,profb,aln->s[a],len_b,len_a,aln->nsip[b]);
84
map[c] = mirror_path(map[c]);
85
}
86
}else{
87
if(aln->nsip[b] == 1){
88
map[c] = ps_dyn(map[c],dp,profa,aln->s[b],len_a,len_b,aln->nsip[a]);
89
}else{
90
if (len_a > len_b){
91
map[c] = pp_dyn(map[c],dp,profa,profb,len_a,len_b);
92
}else{
93
map[c] = pp_dyn(map[c],dp,profb,profa,len_b,len_a);
94
map[c] = mirror_path(map[c]);
95
}
96
}
97
}
98
99
profile[c] = malloc(sizeof(float)*64*(len_a+len_b+2));
100
101
profile[c] = update(profile[a],profile[b],profile[c],map[c],aln->nsip[a],aln->nsip[b]);
102
103
104
aln->sl[c] = map[c][0];
105
106
aln->nsip[c] = aln->nsip[a] + aln->nsip[b];
107
aln->sip[c] = malloc(sizeof(int)*(aln->nsip[a] + aln->nsip[b]));
108
g =0;
109
for (j = aln->nsip[a];j--;){
110
aln->sip[c][g] = aln->sip[a][j];
111
g++;
112
}
113
for (j = aln->nsip[b];j--;){
114
aln->sip[c][g] = aln->sip[b][j];
115
g++;
116
}
117
free(profile[a]);
118
free(profile[b]);
119
}
120
fprintf(stderr,"\r%8.0f percent done\n",100.0);
121
free(profile[numprofiles-1]);
122
free(profile);
123
124
dp_matrix_free(dp);
125
for (i = 32;i--;){
126
free(submatrix[i]);
127
}
128
free(submatrix);
129
return map;
130
}
131
132
/*
133
int** aa_alignment(struct alignment* aln,int* tree,int**submatrix, int** map,int mmbonus)
134
{
135
struct dp_matrix *dp = 0;
136
int i,j,g,a,b,c;
137
int len_a;
138
int len_b;
139
int** profile = 0;
140
int* profa = 0;
141
int* profb = 0;
142
143
int pbonus = 0;
144
145
profile = malloc(sizeof(int*)*numprofiles);
146
for ( i = 0;i< numprofiles;i++){
147
profile[i] = 0;
148
}
149
150
map = malloc(sizeof(int*)*numprofiles);
151
for ( i = 0;i < numprofiles;i++){
152
map[i] = 0;
153
}
154
155
156
dp = dp_matrix_alloc(dp,511,511);
157
c = numseq;
158
for (i = 0; i < (numseq-1);i++){
159
a = tree[i*3];
160
b = tree[i*3+1];
161
c = tree[i*3+2];
162
fprintf(stderr,"Aligning:%d %d->%d\n",a,b,c);
163
len_a = aln->sl[a];
164
len_b = aln->sl[b];
165
dp = dp_matrix_realloc(dp,len_a,len_b);
166
167
map[c] = malloc(sizeof(int) * (len_a+len_b+2));
168
for (j = len_a+len_b+2;j--;){
169
map[c][j] = 0;
170
}
171
if (a < numseq){
172
profile[a] = make_profile(profile[a],aln->s[a],len_a,submatrix);
173
}
174
if (b < numseq){
175
profile[b] = make_profile(profile[b],aln->s[b],len_b,submatrix);
176
}
177
profa = profile[a];
178
profb = profile[b];
179
180
set_gap_penalties(profa,len_a,aln->nsip[b]);
181
set_gap_penalties(profb,len_b,aln->nsip[a]);
182
183
pbonus = mmbonus * aln->nsip[a] * aln->nsip[b];
184
185
if (len_a > len_b){
186
map[c] = aapp_dyn(map[c],dp,profa,profb,len_a,len_b,pbonus);
187
}else{
188
map[c] = aapp_dyn(map[c],dp,profb,profa,len_b,len_a,pbonus);
189
map[c] = mirror_path(map[c]);
190
}
191
192
profile[c] = malloc(sizeof(int)*64*(len_a+len_b+2));
193
profile[c] = update(profa,profb,profile[c],map[c],aln->nsip[a],aln->nsip[b]);
194
195
aln->sl[c] = map[c][0];
196
197
aln->nsip[c] = aln->nsip[a] + aln->nsip[b];
198
aln->sip[c] = malloc(sizeof(int)*(aln->nsip[a] + aln->nsip[b]));
199
g =0;
200
for (j = aln->nsip[a];j--;){
201
aln->sip[c][g] = aln->sip[a][j];
202
g++;
203
}
204
for (j = aln->nsip[b];j--;){
205
aln->sip[c][g] = aln->sip[b][j];
206
g++;
207
}
208
free(profa);
209
free(profb);
210
}
211
212
free(profile[numprofiles-1]);
213
free(profile);
214
215
dp_matrix_free(dp);
216
for (i = 32;i--;){
217
free(submatrix[i]);
218
}
219
free(submatrix);
220
return map;
221
}*/
222
223
/*
224
int** alter_gaps_alignment(struct alignment* aln,int* tree,int**submatrix, int** map,int n,float range,int weight)
225
{
226
struct dp_matrix *dp = 0;
227
int i,j,g,a,b,c;
228
229
int org_gpo = gpo;
230
int org_gpe = gpe;
231
int org_tgpe = tgpe;
232
233
float gpo_step = 0;
234
float gpe_step = 0;
235
float tgpe_step = 0;
236
237
int len_a;
238
int len_b;
239
int** profile = 0;
240
int* profa = 0;
241
int* profb = 0;
242
int* path = 0;
243
244
int* fprofa = 0;
245
int* fprofb = 0;
246
247
if(!(n &1)){
248
n--;
249
}
250
251
float per = 0.0;
252
253
per =(float) range*2/(n+1);
254
255
gpo_step = (float)gpo * per;
256
gpe_step = (float)gpe * per;
257
tgpe_step = (float)tgpe * per;
258
259
260
profile = malloc(sizeof(int*)*numprofiles);
261
for ( i = 0;i< numprofiles;i++){
262
profile[i] = 0;
263
}
264
265
map = malloc(sizeof(int*)*numprofiles);
266
for ( i = 0;i < numprofiles;i++){
267
map[i] = 0;
268
}
269
270
dp = dp_matrix_alloc(dp,511,511);
271
c = numseq;
272
273
for (i = 0; i < (numseq-1);i++){
274
a = tree[i*3];
275
b = tree[i*3+1];
276
c = tree[i*3+2];
277
fprintf(stderr,"Aligning:%d %d->%d\n",a,b,c);
278
len_a = aln->sl[a];
279
len_b = aln->sl[b];
280
dp = dp_matrix_realloc(dp,len_a,len_b);
281
282
map[c] = malloc(sizeof(int) * (len_a+len_b+2));
283
for (j = len_a+len_b+2;j--;){
284
map[c][j] = 0;
285
}
286
if (a < numseq){
287
profile[a] = make_profile(profile[a],aln->s[a],len_a,submatrix);
288
}
289
if (b < numseq){
290
profile[b] = make_profile(profile[b],aln->s[b],len_b,submatrix);
291
}
292
profa = profile[a];
293
profb = profile[b];
294
295
fprofa = malloc(sizeof(int)*(len_a+1)*2);
296
for (j = 0;j < (len_a+1)*2;j++){
297
fprofa[j] = 0;
298
}
299
fprofb = malloc(sizeof(int)*(len_b+1)*2);
300
for (j = 0;j < (len_b+1)*2;j++){
301
fprofb[j] = 0;
302
}
303
304
gpo = org_gpo - ((int)gpo_step* (n/2));
305
gpe = org_gpe - ((int)gpe_step* (n/2));
306
tgpe = org_tgpe - ((int)tgpe_step* (n/2));
307
308
for (j = 0; j < n;j++){
309
set_gap_penalties(profa,len_a,aln->nsip[b]);
310
set_gap_penalties(profb,len_b,aln->nsip[a]);
311
312
path = malloc(sizeof(int) * (len_a+len_b+2));
313
for (g = len_a+len_b+2;g--;){
314
path[g] = 0;
315
}
316
317
if(aln->nsip[a] == 1){
318
if(aln->nsip[b] == 1){
319
path = ss_dyn(submatrix,path,dp,aln->s[a],aln->s[b],len_a,len_b);
320
}else{
321
path = ps_dyn(path,dp,profb,aln->s[a],len_b,len_a,aln->nsip[b]);
322
path = mirror_path(path);
323
}
324
}else{
325
if(aln->nsip[b] == 1){
326
path = ps_dyn(path,dp,profa,aln->s[b],len_a,len_b,aln->nsip[a]);
327
}else{
328
if (len_a > len_b){
329
path = pp_dyn(path,dp,profa,profb,len_a,len_b);
330
}else{
331
path = pp_dyn(path,dp,profb,profa,len_b,len_a);
332
path = mirror_path(path);
333
}
334
}
335
}
336
fprintf(stderr,"Test alignment with gpo:%d gpe:%d tgpe:%d\n",gpo,gpe,tgpe);
337
338
339
add_feature_information_from_alignment(path,fprofa,fprofb,weight/n);
340
341
342
gpo += (int)gpo_step;
343
gpe += (int)gpe_step;
344
tgpe += (int)tgpe_step;
345
}
346
gpo = org_gpo;
347
gpe = org_gpe;
348
tgpe = org_tgpe;
349
350
set_gap_penalties(profa,len_a,aln->nsip[b]);
351
set_gap_penalties(profb,len_b,aln->nsip[a]);
352
353
354
355
if (len_a > len_b){
356
// map[c] = f_only_pp_dyn(map[c],dp,fprofa,fprofb,len_a,len_b,1,2);
357
map[c] = fpp_dyn(map[c],dp,profa,profb,fprofa,fprofb,len_a,len_b,1,2);
358
}else{
359
// map[c] = f_only_pp_dyn(map[c],dp,fprofb,fprofa,len_b,len_a,1,2);
360
map[c] = fpp_dyn(map[c],dp,profb,profa,fprofb,fprofa,len_b,len_a,1,2);
361
map[c] = mirror_path(map[c]);
362
}
363
364
profile[c] = malloc(sizeof(int)*64*(len_a+len_b+2));
365
profile[c] = update(profa,profb,profile[c],map[c],aln->nsip[a],aln->nsip[b]);
366
367
368
aln->sl[c] = map[c][0];
369
370
aln->nsip[c] = aln->nsip[a] + aln->nsip[b];
371
aln->sip[c] = malloc(sizeof(int)*(aln->nsip[a] + aln->nsip[b]));
372
g =0;
373
for (j = aln->nsip[a];j--;){
374
aln->sip[c][g] = aln->sip[a][j];
375
g++;
376
}
377
for (j = aln->nsip[b];j--;){
378
aln->sip[c][g] = aln->sip[b][j];
379
g++;
380
}
381
free(profa);
382
free(profb);
383
384
free(fprofa);
385
free(fprofb);
386
}
387
388
free(profile[numprofiles-1]);
389
free(profile);
390
391
dp_matrix_free(dp);
392
for (i = 32;i--;){
393
free(submatrix[i]);
394
}
395
free(submatrix);
396
return map;
397
}*/
398
399
/*
400
int** test_alignment(struct alignment* aln,int* tree,float **submatrix, int** map,float internal_gap_weight,int window,float strength)
401
{
402
struct dp_matrix *dp = 0;
403
int i,j,g,a,b,c;
404
int len_a;
405
int len_b;
406
float** profile = 0;
407
float* profa = 0;
408
float* profb = 0;
409
410
profile = malloc(sizeof(float*)*numprofiles);
411
for ( i = 0;i< numprofiles;i++){
412
profile[i] = 0;
413
}
414
415
map = malloc(sizeof(int*)*numprofiles);
416
for ( i = 0;i < numprofiles;i++){
417
map[i] = 0;
418
}
419
420
421
dp = dp_matrix_alloc(dp,511,511);
422
c = numseq;
423
for (i = 0; i < (numseq-1);i++){
424
a = tree[i*3];
425
b = tree[i*3+1];
426
c = tree[i*3+2];
427
fprintf(stderr,"Aligning:%d %d->%d\n",a,b,c);
428
len_a = aln->sl[a];
429
len_b = aln->sl[b];
430
dp = dp_matrix_realloc(dp,len_a,len_b);
431
432
map[c] = malloc(sizeof(int) * (len_a+len_b+2));
433
for (j = len_a+len_b+2;j--;){
434
map[c][j] = 0;
435
}
436
if (a < numseq){
437
profile[a] = make_profile2(profile[a],aln->s[a],len_a,submatrix);
438
}
439
if (b < numseq){
440
profile[b] = make_profile2(profile[b],aln->s[b],len_b,submatrix);
441
}
442
profa = profile[a];
443
profb = profile[b];
444
445
set_gap_penalties2(profa,len_a,aln->nsip[b],window,strength);
446
set_gap_penalties2(profb,len_b,aln->nsip[a],window,strength);
447
448
if(aln->nsip[a] == 1){
449
if(aln->nsip[b] == 1){
450
map[c] = ss_dyn2(submatrix,map[c],dp,aln->s[a],aln->s[b],len_a,len_b);
451
}else{
452
// map[c] = ps_dyn2(map[c],dp,profb,aln->s[a],len_b,len_a,aln->nsip[b]);
453
454
map[c] = pp_dyn2(map[c],dp,profb,profa,len_b,len_a);
455
map[c] = mirror_path(map[c]);
456
}
457
}else{
458
if(aln->nsip[b] == 1){
459
// map[c] = ps_dyn2(map[c],dp,profa,aln->s[b],len_a,len_b,aln->nsip[a]);
460
map[c] = pp_dyn2(map[c],dp,profa,profb,len_a,len_b);
461
}else{
462
if (len_a > len_b){
463
map[c] = pp_dyn2(map[c],dp,profa,profb,len_a,len_b);
464
}else{
465
map[c] = pp_dyn2(map[c],dp,profb,profa,len_b,len_a);
466
map[c] = mirror_path(map[c]);
467
}
468
}
469
}
470
471
profile[c] = malloc(sizeof(float)*64*(len_a+len_b+2));
472
profile[c] = update2(profa,profb,profile[c],map[c],aln->nsip[a],aln->nsip[b],internal_gap_weight);
473
474
475
aln->sl[c] = map[c][0];
476
477
aln->nsip[c] = aln->nsip[a] + aln->nsip[b];
478
aln->sip[c] = malloc(sizeof(int)*(aln->nsip[a] + aln->nsip[b]));
479
g =0;
480
for (j = aln->nsip[a];j--;){
481
aln->sip[c][g] = aln->sip[a][j];
482
g++;
483
}
484
for (j = aln->nsip[b];j--;){
485
aln->sip[c][g] = aln->sip[b][j];
486
g++;
487
}
488
free(profa);
489
free(profb);
490
}
491
492
free(profile[numprofiles-1]);
493
free(profile);
494
495
dp_matrix_free(dp);
496
for (i = 32;i--;){
497
free(submatrix[i]);
498
}
499
free(submatrix);
500
return map;
501
}*/
502
/*
503
int** feature_alignment(struct alignment* aln,int* tree,int**submatrix, int** map,struct feature_matrix* fm)
504
{
505
struct dp_matrix *dp = 0;
506
int i,j,g,a,b,c;
507
int len_a;
508
int len_b;
509
int** profile = 0;
510
int* profa = 0;
511
int* profb = 0;
512
513
int** fprofile = 0;
514
int* fprofa = 0;
515
int* fprofb = 0;
516
517
profile = malloc(sizeof(int*)*numprofiles);
518
for ( i = 0;i< numprofiles;i++){
519
profile[i] = 0;
520
}
521
522
fprofile = malloc(sizeof(int*)*numprofiles);
523
for ( i = 0;i< numprofiles;i++){
524
fprofile[i] = 0;
525
}
526
527
map = malloc(sizeof(int*)*numprofiles);
528
for ( i = 0;i < numprofiles;i++){
529
map[i] = 0;
530
}
531
532
dp = dp_matrix_alloc(dp,511,511);
533
c = numseq;
534
//if(!param->dna){
535
for (i = 0; i < (numseq-1);i++){
536
a = tree[i*3];
537
b = tree[i*3+1];
538
c = tree[i*3+2];
539
fprintf(stderr,"Aligning:%d %d->%d\n",a,b,c);
540
len_a = aln->sl[a];
541
len_b = aln->sl[b];
542
dp = dp_matrix_realloc(dp,len_a,len_b);
543
544
map[c] = malloc(sizeof(int) * (len_a+len_b+2));
545
for (j = len_a+len_b+2;j--;){
546
map[c][j] = 0;
547
}
548
if (a < numseq){
549
profile[a] = make_profile(profile[a],aln->s[a],len_a,submatrix);
550
// fprintf(stderr,"Making feature profile for %d (%s)\n",a,aln->sn[a]);
551
fprofile[a] = make_feature_profile(fprofile[a],aln->ft[a],len_a,fm);
552
}
553
if (b < numseq){
554
profile[b] = make_profile(profile[b],aln->s[b],len_b,submatrix);
555
// fprintf(stderr,"Making feature profile for %d (%s)\n",b,aln->sn[b]);
556
fprofile[b] = make_feature_profile(fprofile[b],aln->ft[b],len_b,fm);
557
}
558
//profa = profile[a];
559
//profb = profile[b];
560
profa = profile[a]+64;
561
profb = profile[b]+64;
562
563
fprofa = fprofile[a];
564
fprofb = fprofile[b];
565
566
set_gap_penalties(profile[a],len_a,aln->nsip[b]);
567
set_gap_penalties(profile[b],len_b,aln->nsip[a]);
568
569
if (len_a > len_b){
570
map[c] = fpp_dyn(map[c],dp,profa,profb,fprofa,fprofb,len_a,len_b,fm->mdim,fm->stride);
571
}else{
572
map[c] = fpp_dyn(map[c],dp,profb,profa,fprofb,fprofa,len_b,len_a,fm->mdim,fm->stride);
573
map[c] = mirror_path(map[c]);
574
}
575
576
profile[c] = malloc(sizeof(int)*64*(len_a+len_b+2));
577
profile[c] = update(profile[a],profile[b],profile[c],map[c],aln->nsip[a],aln->nsip[b]);
578
579
fprofile[c] = malloc(sizeof(int)*fm->stride*(len_a+len_b+2));
580
fprofile[c] = feature_update(fprofa,fprofb,fprofile[c],map[c],fm->stride);
581
582
aln->sl[c] = map[c][0];
583
584
aln->nsip[c] = aln->nsip[a] + aln->nsip[b];
585
aln->sip[c] = malloc(sizeof(int)*(aln->nsip[a] + aln->nsip[b]));
586
g =0;
587
for (j = aln->nsip[a];j--;){
588
aln->sip[c][g] = aln->sip[a][j];
589
g++;
590
}
591
for (j = aln->nsip[b];j--;){
592
aln->sip[c][g] = aln->sip[b][j];
593
g++;
594
}
595
free(profile[a]);
596
free(profile[b]);
597
598
free(fprofa);
599
free(fprofb);
600
601
}
602
603
free(profile[numprofiles-1]);
604
free(profile);
605
606
free(fprofile[numprofiles-1]);
607
free(fprofile );
608
609
dp_matrix_free(dp);
610
for (i = 32;i--;){
611
free(submatrix[i]);
612
}
613
free(submatrix);
614
free_feature_matrix(fm);
615
return map;
616
}*/
617
618
struct ntree_data* ntree_sub_alignment(struct ntree_data* ntree_data,int* tree,int num)
619
{
620
struct dp_matrix *dp = 0;
621
struct alignment* aln = 0;
622
int i,j,g,a,b,c;
623
int len_a;
624
int len_b;
625
float** local_profile = 0;
626
float* profa = 0;
627
float* profb = 0;
628
629
int** local_map = 0;
630
int* local_sl = 0;
631
int* local_nsip = 0;
632
int** local_sip = 0;
633
634
int* which_to_alloc = 0;
635
636
aln = ntree_data->aln;
637
638
which_to_alloc = malloc(sizeof(int*)*numprofiles);
639
for ( i = 0;i< numprofiles;i++){
640
which_to_alloc[i] = 0;
641
}
642
643
local_profile = malloc(sizeof(float*)*numprofiles);
644
local_sl = malloc(sizeof(int)*numprofiles);
645
local_nsip = malloc(sizeof(int)*numprofiles);
646
local_sip = malloc(sizeof(int*)*numprofiles);
647
648
649
for (i = 0; i < num-1;i++){
650
a = tree[i*3+1];
651
if(!which_to_alloc[a]){
652
which_to_alloc[a] = 1;
653
}
654
b = tree[i*3+2];
655
if(!which_to_alloc[b]){
656
which_to_alloc[b] = 1;
657
}
658
c = tree[i*3+3];
659
if(!which_to_alloc[c]){
660
which_to_alloc[c] = 2;
661
}
662
}
663
//for ( i = 0;i< numprofiles;i++){
664
// fprintf(stderr,"alloc?:%d %d\n",i,which_to_alloc[i]);
665
//}
666
667
// exit(0);
668
for ( i = 0;i< numprofiles;i++){
669
if(which_to_alloc[i] == 1){
670
local_profile[i] = ntree_data->profile[i];
671
local_sl[i] = aln->sl[i];
672
local_nsip[i] = aln->nsip[i];
673
local_sip[i] = malloc(sizeof(int*)*aln->nsip[i]);
674
for(j = 0;j < aln->nsip[i];j++){
675
local_sip[i][j] = aln->sip[i][j];
676
}
677
}else{
678
local_profile[i] = 0;
679
local_sl[i] = 0;
680
local_nsip[i] = 0;
681
local_sip[i] = 0;
682
}
683
}
684
/*
685
for ( i = 0;i< numprofiles;i++){
686
local_profile[i] = ntree_data->profile[i];
687
local_sl[i] = aln->sl[i];
688
local_nsip[i] = aln->nsip[i];
689
if(aln->sip[i]){
690
fprintf(stderr,"Allocing..:%d\n",aln->nsip[i]);
691
local_sip[i] = malloc(sizeof(int*)*aln->nsip[i]);
692
for(j = 0;j < aln->nsip[i];j++){
693
local_sip[i][j] = aln->sip[i][j];
694
}
695
}else{
696
local_sip[i] = 0;
697
}
698
}*/
699
700
local_map = malloc(sizeof(int*)*numprofiles);
701
for ( i = 0;i < numprofiles;i++){
702
local_map[i] = 0;
703
}
704
705
706
dp = dp_matrix_alloc(dp,511,511);
707
c = numseq;
708
for (i = 0; i < num-1;i++){
709
a = tree[i*3+1];
710
b = tree[i*3+2];
711
c = tree[i*3+3];
712
// fprintf(stderr,"Aligning:%d %d->%d\n",a,b,c);
713
len_a = local_sl[a];
714
len_b = local_sl[b];
715
dp = dp_matrix_realloc(dp,len_a,len_b);
716
717
local_map[c] = malloc(sizeof(int) * (len_a+len_b+2));
718
for (j = len_a+len_b+2;j--;){
719
local_map[c][j] = 0;
720
}
721
if (a < numseq){
722
local_profile[a] = make_profile(local_profile[a],aln->s[a],len_a,ntree_data->submatrix);
723
}
724
if (b < numseq){
725
local_profile[b] = make_profile(local_profile[b],aln->s[b],len_b,ntree_data->submatrix);
726
}
727
profa = local_profile[a];
728
profb = local_profile[b];
729
730
set_gap_penalties(profa,len_a,local_nsip[b]);
731
set_gap_penalties(profb,len_b,local_nsip[a]);
732
733
if(local_nsip[a] == 1){
734
if(local_nsip[b] == 1){
735
local_map[c] = ss_dyn(ntree_data->submatrix,local_map[c],dp,aln->s[a],aln->s[b],len_a,len_b);
736
}else{
737
local_map[c] = ps_dyn(local_map[c],dp,profb,aln->s[a],len_b,len_a,local_nsip[b]);
738
local_map[c] = mirror_path(local_map[c]);
739
}
740
}else{
741
if(local_nsip[b] == 1){
742
local_map[c] = ps_dyn(local_map[c],dp,profa,aln->s[b],len_a,len_b,local_nsip[a]);
743
}else{
744
if (len_a > len_b){
745
local_map[c] = pp_dyn(local_map[c],dp,profa,profb,len_a,len_b);
746
}else{
747
local_map[c] = pp_dyn(local_map[c],dp,profb,profa,len_b,len_a);
748
local_map[c] = mirror_path(local_map[c]);
749
}
750
}
751
}
752
753
local_profile[c] = malloc(sizeof(float)*64*(len_a+len_b+2));
754
local_profile[c] = update(profa,profb,local_profile[c],local_map[c],local_nsip[a],local_nsip[b]);
755
756
local_sl[c] = local_map[c][0];
757
758
local_nsip[c] = local_nsip[a] + local_nsip[b];
759
local_sip[c] = malloc(sizeof(int)*(local_nsip[a] + local_nsip[b]));
760
g =0;
761
for (j = local_nsip[a];j--;){
762
local_sip[c][g] = local_sip[a][j];
763
g++;
764
}
765
for (j = local_nsip[b];j--;){
766
local_sip[c][g] = local_sip[b][j];
767
g++;
768
}
769
// free(profa);
770
// free(profb);
771
}
772
773
if(ntree_data->profile[c]){
774
if(ntree_data->map[c][ntree_data->map[c][0]+2] < local_map[c][local_map[c][0]+2]){
775
fprintf(stderr,"%d\n",local_map[c][local_map[c][0]+2]);
776
//remove old map,profile,etc..
777
for (i = 0; i < num-1;i++){
778
c = tree[i*3+3];
779
free(ntree_data->map[c]);
780
free(ntree_data->profile[c]);
781
free(aln->sip[c]);
782
ntree_data->map[c] = malloc(sizeof(int)*(local_map[c][0]+3));
783
for (j = 0; j < local_map[c][0]+3;j++){
784
ntree_data->map[c][j] = local_map[c][j];
785
}
786
aln->sip[c] = malloc(sizeof(int)*local_nsip[c]);
787
aln->nsip[c] = local_nsip[c];
788
for (j = 0; j < local_nsip[c];j++){
789
aln->sip[c][j] = local_sip[c][j];
790
}
791
aln->sl[c] = local_sl[c];
792
793
}
794
ntree_data->profile[c] = malloc(sizeof(int)*64*(aln->sl[c]+1));
795
for (i = 0; i < (64*(aln->sl[c]+1));i++){
796
ntree_data->profile[c][i] = local_profile[c][i];
797
}
798
ntree_data->tree[0] -= (tree[0]-1);
799
for (j = 1; j < tree[0];j++){
800
ntree_data->tree[ntree_data->tree[0]+j-1] = tree[j];
801
}
802
ntree_data->tree[0] += (tree[0]-1);
803
804
}else{
805
fprintf(stderr,"no improvement\n");
806
}
807
}else{
808
fprintf(stderr,"%d\n",local_map[c][local_map[c][0]+2]);
809
for (i = 0; i < num-1;i++){
810
c = tree[i*3+3];
811
ntree_data->map[c] = malloc(sizeof(int)*(local_map[c][0]+3));
812
for (j = 0; j < local_map[c][0]+3;j++){
813
ntree_data->map[c][j] = local_map[c][j];
814
}
815
816
aln->sip[c] = malloc(sizeof(int)*local_nsip[c]);
817
aln->nsip[c] = local_nsip[c];
818
for (j = 0; j < local_nsip[c];j++){
819
aln->sip[c][j] = local_sip[c][j];
820
}
821
aln->sl[c] = local_sl[c];
822
}
823
ntree_data->profile[c] = malloc(sizeof(int)*64*(aln->sl[c]+1));
824
for (i = 0; i < (64*(aln->sl[c]+1));i++){
825
ntree_data->profile[c][i] = local_profile[c][i];
826
}
827
for (j = 1; j < tree[0];j++){
828
ntree_data->tree[ntree_data->tree[0]+j-1] = tree[j];
829
}
830
ntree_data->tree[0] += tree[0]-1;
831
}
832
833
for ( i = 0;i< numprofiles;i++){
834
if(which_to_alloc[i] == 1){
835
free(local_sip[i]);
836
if(i < numseq){
837
free(local_profile[i]);
838
}
839
}
840
if(which_to_alloc[i] == 2){
841
free(local_profile[i]);
842
free(local_map[i]);
843
free(local_sip[i]);
844
}
845
846
}
847
848
free(which_to_alloc);
849
free(local_map);
850
free(local_sip);
851
free(local_nsip);
852
free(local_profile);
853
free(local_sl);
854
855
dp_matrix_free(dp);
856
return ntree_data;
857
}
858
859
struct ntree_data* ntree_alignment(struct ntree_data* ntree_data)
860
{
861
int i;
862
ntree_data->profile = malloc(sizeof(float*)*numprofiles);
863
for ( i = 0;i< numprofiles;i++){
864
ntree_data->profile[i] = 0;
865
}
866
867
ntree_data->map = malloc(sizeof(int*)*numprofiles);
868
for ( i = 0;i < numprofiles;i++){
869
ntree_data->map[i] = 0;
870
}
871
872
ntree_data = alignntree(ntree_data,ntree_data->realtree);
873
874
for ( i = 0;i< numprofiles;i++){
875
if(ntree_data->profile[i]){
876
free(ntree_data->profile[i]);
877
}
878
}
879
free(ntree_data->profile);
880
881
for (i = 32;i--;){
882
free(ntree_data->submatrix[i]);
883
}
884
free(ntree_data->submatrix);
885
free_real_tree(ntree_data->realtree);
886
return ntree_data;
887
}
888
Loggerhead is a web-based interface for Breezy
Version: 2.0.1