~ubuntu-branches/debian/jessie/arb/jessie

« back to all changes in this revision

Viewing changes to GDE/MAFFT/mafft-7.055-with-extensions/core/genGalign11.c

  • Committer: Package Import Robot
  • Author(s): Elmar Pruesse, Andreas Tille, Elmar Pruesse
  • Date: 2014-09-02 15:15:06 UTC
  • mfrom: (1.1.6)
  • Revision ID: package-import@ubuntu.com-20140902151506-jihq58b3iz342wif
Tags: 6.0.2-1
[ Andreas Tille ]
* New upstream version
  Closes: #741890
* debian/upstream -> debian/upstream/metadata
* debian/control:
   - Build-Depends: added libglib2.0-dev
   - Depends: added mafft, mrbayes
* debian/rules
   - Add explicite --remove-section=.comment option to manual strip call
* cme fix dpkg-control
* arb-common.dirs: Do not create unneeded lintian dir
* Add turkish debconf translation (thanks for the patch to Mert Dirik
  <mertdirik@gmail.com>)
  Closes: #757497

[ Elmar Pruesse ]
* patches removed:
   - 10_config.makefiles.patch,
     80_no_GL.patch
       removed in favor of creating file from config.makefile.template via 
       sed in debian/control
   - 20_Makefile_main.patch
       merged upstream
   - 21_Makefiles.patch
       no longer needed
   - 30_tmpfile_CVE-2008-5378.patch: 
       merged upstream
   - 50_fix_gcc-4.8.patch:
       merged upstream
   - 40_add_libGLU.patch:
       libGLU not needed for arb_ntree)
   - 60_use_debian_packaged_raxml.patch:
       merged upstream
   - 70_hardening.patch
       merged upstream
   - 72_add_math_lib_to_linker.patch
       does not appear to be needed
* patches added:
   - 10_upstream_r12793__show_db_load_progress:
       backported patch showing progress while ARB is loading a database
       (needed as indicator/splash screen while ARB is launching)
   - 20_upstream_r12794__socket_permissions:
       backported security fix
   - 30_upstream_r12814__desktop_keywords:
       backported add keywords to desktop (fixes lintian warning)
   - 40_upstream_r12815__lintian_spelling:
       backported fix for lintian reported spelling errors
   - 50_private_nameservers
       change configuration to put nameservers into users home dirs
       (avoids need for shared writeable directory)
   - 60_use_debian_phyml
       use phyml from debian package for both interfaces in ARB
* debian/rules:
   - create config.makefile from override_dh_configure target
   - use "make tarfile" in override_dh_install
   - remove extra cleaning not needed for ARB 6
   - use "dh_install --list-missing" to avoid missing files
   - added override_dh_fixperms target
* debian/control:
   - added libarb-dev package
   - Depends: added phyml, xdg-utils
   - Suggests: removed phyml
   - fix lintian duplicate-short-description (new descriptions)
* debian/*.install:
   - "unrolled" confusing globbing to select files
   - pick files from debian/tmp
   - moved all config files to /etc/arb
* debian/arb-common.templates: updated
* scripts:
   - removed arb-add-pt-server
   - launch-wrapper: 
     - only add demo.arb to newly created $ARBUSERDATA
     - pass commandline arguments through bin/arb wrapper
   - preinst: removing old PT server index files on upgrade from 5.5*
   - postinst: set setgid on shared PT dir
* rewrote arb.1 manfile
* added file icon for ARB databases
* using upstream arb_tcp.dat

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include "mltaln.h"
 
2
#include "dp.h"
 
3
 
 
4
#define DEBUG 0
 
5
#define XXXXXXX    0
 
6
#define USE_PENALTY_EX  1
 
7
 
 
8
 
 
9
#if 1
 
10
static void match_calc( float *match, char **s1, char **s2, int i1, int lgth2 ) 
 
11
{
 
12
        char *seq2 = s2[0];
 
13
        int *intptr = amino_dis[(int)s1[0][i1]];
 
14
 
 
15
        while( lgth2-- )
 
16
                *match++ = intptr[(int)*seq2++];
 
17
}
 
18
#else
 
19
static void match_calc( float *match, char **s1, char **s2, int i1, int lgth2 )
 
20
{
 
21
        int j;
 
22
 
 
23
        for( j=0; j<lgth2; j++ )
 
24
                match[j] = amino_dis[(*s1)[i1]][(*s2)[j]];
 
25
}
 
26
#endif
 
27
 
 
28
static float genGtracking( float *lasthorizontalw, float *lastverticalw, 
 
29
                                                char **seq1, char **seq2, 
 
30
                        char **mseq1, char **mseq2, 
 
31
                        float **cpmx1, float **cpmx2, 
 
32
                        int **ijpi, int **ijpj )
 
33
{
 
34
        int i, j, l, iin, jin, ifi, jfi, lgth1, lgth2, k, limk;
 
35
//      char gap[] = "-";
 
36
        char *gap;
 
37
        gap = newgapstr;
 
38
        lgth1 = strlen( seq1[0] );
 
39
        lgth2 = strlen( seq2[0] );
 
40
 
 
41
 
 
42
#if 0
 
43
        for( i=0; i<lgth1; i++ ) 
 
44
        {
 
45
                fprintf( stderr, "lastverticalw[%d] = %f\n", i, lastverticalw[i] );
 
46
        }
 
47
#endif
 
48
 
 
49
    for( i=0; i<lgth1+1; i++ ) 
 
50
    {
 
51
        ijpi[i][0] = -1;
 
52
                ijpj[i][0] = -1; // ???
 
53
    }
 
54
    for( j=0; j<lgth2+1; j++ ) 
 
55
    {
 
56
                ijpi[0][j] = -1; // ???
 
57
        ijpj[0][j] = -1;
 
58
    }
 
59
 
 
60
 
 
61
        mseq1[0] += lgth1+lgth2;
 
62
        *mseq1[0] = 0;
 
63
        mseq2[0] += lgth1+lgth2;
 
64
        *mseq2[0] = 0;
 
65
        iin = lgth1; jin = lgth2;
 
66
        limk = lgth1+lgth2 + 1;
 
67
        for( k=0; k<limk; k++ ) 
 
68
        {
 
69
                ifi = ( ijpi[iin][jin] );
 
70
                jfi = ( ijpj[iin][jin] );
 
71
                l = iin - ifi;
 
72
                while( --l ) 
 
73
                {
 
74
                        *--mseq1[0] = seq1[0][ifi+l];
 
75
                        *--mseq2[0] = *gap;
 
76
                        k++;
 
77
                }
 
78
                l= jin - jfi;
 
79
                while( --l )
 
80
                {
 
81
                        *--mseq1[0] = *gap;
 
82
                        *--mseq2[0] = seq2[0][jfi+l];
 
83
                        k++;
 
84
                }
 
85
                if( iin <= 0 || jin <= 0 ) break;
 
86
                *--mseq1[0] = seq1[0][ifi];
 
87
                *--mseq2[0] = seq2[0][jfi];
 
88
 
 
89
//              fprintf( stdout, "mseq1 = %s\n", mseq1[0] );
 
90
//              fprintf( stdout, "mseq2 = %s\n", mseq2[0] );
 
91
 
 
92
                k++;
 
93
                iin = ifi; jin = jfi;
 
94
        }
 
95
        return( 0.0 );
 
96
}
 
97
 
 
98
 
 
99
float genG__align11( char **seq1, char **seq2, int alloclen )
 
100
/* score no keisan no sai motokaraaru gap no atukai ni mondai ga aru */
 
101
{
 
102
//      int k;
 
103
        register int i, j;
 
104
        int lasti;                      /* outgap == 0 -> lgth1, outgap == 1 -> lgth1+1 */
 
105
        int lgth1, lgth2;
 
106
        int resultlen;
 
107
        float wm;   /* int ?????? */
 
108
        float g;
 
109
        float *currentw, *previousw;
 
110
        float fpenalty = (float)penalty;
 
111
        float fpenalty_OP = (float)penalty_OP;
 
112
#if USE_PENALTY_EX
 
113
        float fpenalty_ex = (float)penalty_ex;
 
114
#endif
 
115
#if 1
 
116
        float *wtmp;
 
117
        int *ijpipt;
 
118
        int *ijpjpt;
 
119
        float *mjpt, *Mjpt, *prept, *curpt;
 
120
        int *mpjpt, *Mpjpt;
 
121
#endif
 
122
        static float mi, *m;
 
123
        static float Mi, *largeM;
 
124
        static int **ijpi;
 
125
        static int **ijpj;
 
126
        static int mpi, *mp;
 
127
        static int Mpi, *Mp;
 
128
        static float *w1, *w2;
 
129
        static float *match;
 
130
        static float *initverticalw;    /* kufuu sureba iranai */
 
131
        static float *lastverticalw;    /* kufuu sureba iranai */
 
132
        static char **mseq1;
 
133
        static char **mseq2;
 
134
        static char **mseq;
 
135
        static float **cpmx1;
 
136
        static float **cpmx2;
 
137
        static int **intwork;
 
138
        static float **floatwork;
 
139
        static int orlgth1 = 0, orlgth2 = 0;
 
140
        float tbk;
 
141
        int tbki, tbkj;
 
142
 
 
143
        wm = 0.0;
 
144
 
 
145
        if( orlgth1 == 0 )
 
146
        {
 
147
                mseq1 = AllocateCharMtx( njob, 0 );
 
148
                mseq2 = AllocateCharMtx( njob, 0 );
 
149
        }
 
150
 
 
151
 
 
152
        lgth1 = strlen( seq1[0] );
 
153
        lgth2 = strlen( seq2[0] );
 
154
 
 
155
 
 
156
 
 
157
        if( lgth1 <= 0 || lgth2 <= 0 )
 
158
        {
 
159
                fprintf( stderr, "WARNING (g11): lgth1=%d, lgth2=%d\n", lgth1, lgth2 );
 
160
        }
 
161
 
 
162
        if( lgth1 > orlgth1 || lgth2 > orlgth2 )
 
163
        {
 
164
                int ll1, ll2;
 
165
 
 
166
                if( orlgth1 > 0 && orlgth2 > 0 )
 
167
                {
 
168
                        FreeFloatVec( w1 );
 
169
                        FreeFloatVec( w2 );
 
170
                        FreeFloatVec( match );
 
171
                        FreeFloatVec( initverticalw );
 
172
                        FreeFloatVec( lastverticalw );
 
173
 
 
174
                        FreeFloatVec( m );
 
175
                        FreeIntVec( mp );
 
176
                        FreeFloatVec( largeM );
 
177
                        FreeIntVec( Mp );
 
178
 
 
179
                        FreeCharMtx( mseq );
 
180
 
 
181
 
 
182
                        FreeFloatMtx( cpmx1 );
 
183
                        FreeFloatMtx( cpmx2 );
 
184
 
 
185
                        FreeFloatMtx( floatwork );
 
186
                        FreeIntMtx( intwork );
 
187
                }
 
188
 
 
189
                ll1 = MAX( (int)(1.3*lgth1), orlgth1 ) + 100;
 
190
                ll2 = MAX( (int)(1.3*lgth2), orlgth2 ) + 100;
 
191
 
 
192
#if DEBUG
 
193
                fprintf( stderr, "\ntrying to allocate (%d+%d)xn matrices ... ", ll1, ll2 );
 
194
#endif
 
195
 
 
196
                w1 = AllocateFloatVec( ll2+2 );
 
197
                w2 = AllocateFloatVec( ll2+2 );
 
198
                match = AllocateFloatVec( ll2+2 );
 
199
 
 
200
                initverticalw = AllocateFloatVec( ll1+2 );
 
201
                lastverticalw = AllocateFloatVec( ll1+2 );
 
202
 
 
203
                m = AllocateFloatVec( ll2+2 );
 
204
                mp = AllocateIntVec( ll2+2 );
 
205
                largeM = AllocateFloatVec( ll2+2 );
 
206
                Mp = AllocateIntVec( ll2+2 );
 
207
 
 
208
                mseq = AllocateCharMtx( njob, ll1+ll2 );
 
209
 
 
210
                cpmx1 = AllocateFloatMtx( 26, ll1+2 );
 
211
                cpmx2 = AllocateFloatMtx( 26, ll2+2 );
 
212
 
 
213
                floatwork = AllocateFloatMtx( 26, MAX( ll1, ll2 )+2 ); 
 
214
                intwork = AllocateIntMtx( 26, MAX( ll1, ll2 )+2 ); 
 
215
 
 
216
#if DEBUG
 
217
                fprintf( stderr, "succeeded\n" );
 
218
#endif
 
219
 
 
220
                orlgth1 = ll1 - 100;
 
221
                orlgth2 = ll2 - 100;
 
222
        }
 
223
 
 
224
 
 
225
        mseq1[0] = mseq[0];
 
226
        mseq2[0] = mseq[1];
 
227
 
 
228
 
 
229
        if( orlgth1 > commonAlloc1 || orlgth2 > commonAlloc2 )
 
230
        {
 
231
                int ll1, ll2;
 
232
 
 
233
                if( commonAlloc1 && commonAlloc2 )
 
234
                {
 
235
                        FreeIntMtx( commonIP );
 
236
                        FreeIntMtx( commonJP );
 
237
                }
 
238
 
 
239
                ll1 = MAX( orlgth1, commonAlloc1 );
 
240
                ll2 = MAX( orlgth2, commonAlloc2 );
 
241
 
 
242
#if DEBUG
 
243
                fprintf( stderr, "\n\ntrying to allocate %dx%d matrices ... ", ll1+1, ll2+1 );
 
244
#endif
 
245
 
 
246
                commonIP = AllocateIntMtx( ll1+10, ll2+10 );
 
247
                commonJP = AllocateIntMtx( ll1+10, ll2+10 );
 
248
 
 
249
#if DEBUG
 
250
                fprintf( stderr, "succeeded\n\n" );
 
251
#endif
 
252
 
 
253
                commonAlloc1 = ll1;
 
254
                commonAlloc2 = ll2;
 
255
        }
 
256
        ijpi = commonIP;
 
257
        ijpj = commonJP;
 
258
 
 
259
 
 
260
#if 0
 
261
        for( i=0; i<lgth1; i++ ) 
 
262
                fprintf( stderr, "ogcp1[%d]=%f\n", i, ogcp1[i] );
 
263
#endif
 
264
 
 
265
        currentw = w1;
 
266
        previousw = w2;
 
267
 
 
268
 
 
269
        match_calc( initverticalw, seq2, seq1, 0, lgth1 );
 
270
 
 
271
 
 
272
        match_calc( currentw, seq1, seq2, 0, lgth2 );
 
273
 
 
274
        if( outgap == 1 )
 
275
        {
 
276
                for( i=1; i<lgth1+1; i++ )
 
277
                {
 
278
                        initverticalw[i] += fpenalty;
 
279
                }
 
280
                for( j=1; j<lgth2+1; j++ )
 
281
                {
 
282
                        currentw[j] += fpenalty;
 
283
                }
 
284
        }
 
285
 
 
286
        for( j=1; j<lgth2+1; ++j ) 
 
287
        {
 
288
                m[j] = currentw[j-1]; mp[j] = 0;
 
289
                largeM[j] = currentw[j-1]; Mp[j] = 0;
 
290
        }
 
291
 
 
292
        if( lgth2 == 0 )
 
293
                lastverticalw[0] = 0.0;               // lgth2==0 no toki error
 
294
        else
 
295
                lastverticalw[0] = currentw[lgth2-1]; // lgth2==0 no toki error
 
296
 
 
297
        if( outgap ) lasti = lgth1+1; else lasti = lgth1;
 
298
 
 
299
#if XXXXXXX
 
300
fprintf( stderr, "currentw = \n" );
 
301
for( i=0; i<lgth1+1; i++ )
 
302
{
 
303
        fprintf( stderr, "%5.2f ", currentw[i] );
 
304
}
 
305
fprintf( stderr, "\n" );
 
306
fprintf( stderr, "initverticalw = \n" );
 
307
for( i=0; i<lgth2+1; i++ )
 
308
{
 
309
        fprintf( stderr, "%5.2f ", initverticalw[i] );
 
310
}
 
311
fprintf( stderr, "\n" );
 
312
#endif
 
313
 
 
314
        for( i=1; i<lasti; i++ )
 
315
        {
 
316
                wtmp = previousw; 
 
317
                previousw = currentw;
 
318
                currentw = wtmp;
 
319
 
 
320
                previousw[0] = initverticalw[i-1];
 
321
 
 
322
                match_calc( currentw, seq1, seq2, i, lgth2 );
 
323
#if XXXXXXX
 
324
fprintf( stderr, "\n" );
 
325
fprintf( stderr, "i=%d\n", i );
 
326
fprintf( stderr, "currentw = \n" );
 
327
for( j=0; j<lgth2; j++ )
 
328
{
 
329
        fprintf( stderr, "%5.2f ", currentw[j] );
 
330
}
 
331
fprintf( stderr, "\n" );
 
332
#endif
 
333
#if XXXXXXX
 
334
fprintf( stderr, "\n" );
 
335
fprintf( stderr, "i=%d\n", i );
 
336
fprintf( stderr, "currentw = \n" );
 
337
for( j=0; j<lgth2; j++ )
 
338
{
 
339
        fprintf( stderr, "%5.2f ", currentw[j] );
 
340
}
 
341
fprintf( stderr, "\n" );
 
342
#endif
 
343
                currentw[0] = initverticalw[i];
 
344
 
 
345
                mi = previousw[0]; mpi = 0;
 
346
                Mi = previousw[0]; Mpi = 0;
 
347
 
 
348
                ijpipt = ijpi[i] + 1;
 
349
                ijpjpt = ijpj[i] + 1;
 
350
                mjpt = m + 1;
 
351
                Mjpt = largeM + 1;
 
352
                prept = previousw;
 
353
                curpt = currentw + 1;
 
354
                mpjpt = mp + 1;
 
355
                Mpjpt = Mp + 1;
 
356
                tbk = -9999999.9;
 
357
                tbki = 0;
 
358
                tbkj = 0;
 
359
                for( j=1; j<lgth2+1; j++ )
 
360
                {
 
361
                        wm = *prept;
 
362
                        *ijpipt = i-1;
 
363
                        *ijpjpt = j-1;
 
364
 
 
365
#if 0
 
366
                        fprintf( stderr, "%5.0f->", wm );
 
367
#endif
 
368
#if 0
 
369
                        fprintf( stderr, "%5.0f?", g );
 
370
#endif
 
371
                        if( (g=mi+fpenalty) > wm )
 
372
                        {
 
373
                                wm = g;
 
374
//                              *ijpipt = i - 1; // iranai
 
375
                                *ijpjpt = mpi;
 
376
                        }
 
377
                        if( (g=*prept) >= mi )
 
378
                        {
 
379
                                mi = g;
 
380
                                mpi = j-1;
 
381
                        }
 
382
#if USE_PENALTY_EX
 
383
                        mi += fpenalty_ex;
 
384
#endif
 
385
 
 
386
#if 0 
 
387
                        fprintf( stderr, "%5.0f?", g );
 
388
#endif
 
389
                        if( (g=*mjpt + fpenalty) > wm )
 
390
                        {
 
391
                                wm = g;
 
392
                                *ijpipt = *mpjpt;
 
393
                                *ijpjpt = j - 1; //IRU!
 
394
                        }
 
395
                        if( (g=*prept) >= *mjpt )
 
396
                        {
 
397
                                *mjpt = g;
 
398
                                *mpjpt = i-1;
 
399
                        }
 
400
#if USE_PENALTY_EX
 
401
                        m[j] += fpenalty_ex;
 
402
#endif
 
403
 
 
404
#if 1
 
405
                        g =  tbk + fpenalty_OP;
 
406
                        if( g > wm )
 
407
                        {
 
408
                                wm = g;
 
409
                                *ijpipt = tbki;
 
410
                                *ijpjpt = tbkj;
 
411
//                              fprintf( stderr, "hit! i%d, j%d, ijpi = %d, ijpj = %d\n", i, j, *ijpipt, *ijpjpt );
 
412
                        }
 
413
                        if( Mi > tbk )
 
414
                        {
 
415
                                tbk = Mi; //error desu.
 
416
                                tbki = i-1;
 
417
                                tbkj = Mpi;
 
418
                        }
 
419
                        if( *Mjpt > tbk )
 
420
                        {
 
421
                                tbk = *Mjpt;
 
422
                                tbki = *Mpjpt;
 
423
                                tbkj = j-1;
 
424
                        }
 
425
 
 
426
                        if( *prept > *Mjpt )
 
427
                        {
 
428
                                *Mjpt = *prept;
 
429
                                *Mpjpt = i-1;
 
430
                        }
 
431
                        if( *prept > Mi )
 
432
                        {
 
433
                                Mi = *prept;
 
434
                                Mpi = j-1;
 
435
                        }
 
436
 
 
437
#endif
 
438
 
 
439
 
 
440
#if 0
 
441
                        fprintf( stderr, "%5.0f ", wm );
 
442
#endif
 
443
                        *curpt++ += wm;
 
444
                        ijpipt++;
 
445
                        ijpjpt++;
 
446
                        mjpt++;
 
447
                        Mjpt++;
 
448
                        prept++;
 
449
                        mpjpt++;
 
450
                        Mpjpt++;
 
451
                }
 
452
                lastverticalw[i] = currentw[lgth2-1]; // lgth2==0 no toki error
 
453
        }
 
454
#if 0
 
455
        for( i=0; i<lgth1; i++ )
 
456
        {
 
457
                for( j=0; j<lgth2; j++ )
 
458
                {
 
459
                        fprintf( stdout, "i,j=%d,%d - ijpi,ijpj=%d,%d\n", i, j, ijpi[i][j], ijpj[i][j] );
 
460
                }
 
461
        }
 
462
        fflush( stdout );
 
463
#endif
 
464
 
 
465
        genGtracking( currentw, lastverticalw, seq1, seq2, mseq1, mseq2, cpmx1, cpmx2, ijpi, ijpj );
 
466
 
 
467
        resultlen = strlen( mseq1[0] );
 
468
        if( alloclen < resultlen || resultlen > N )
 
469
        {
 
470
                fprintf( stderr, "alloclen=%d, resultlen=%d, N=%d\n", alloclen, resultlen, N );
 
471
                ErrorExit( "LENGTH OVER!\n" );
 
472
        }
 
473
 
 
474
 
 
475
        strcpy( seq1[0], mseq1[0] );
 
476
        strcpy( seq2[0], mseq2[0] );
 
477
#if 0
 
478
        fprintf( stderr, "\n" );
 
479
        fprintf( stderr, ">\n%s\n", mseq1[0] );
 
480
        fprintf( stderr, ">\n%s\n", mseq2[0] );
 
481
        fprintf( stderr, "wm = %f\n", wm );
 
482
#endif
 
483
 
 
484
        return( wm );
 
485
}
 
486