| /************************/ |
| /* sim.c, 27 Oct 1994 */ |
| /* runtime --> sim.res */ |
| /************************/ |
| |
| /* A PROGRAM FOR LOCAL SIMILARITIES WITH AFFINE WEIGHTS: |
| |
| Modified 12/6/90 for consistency with the "blast", "convert" and "lav" |
| programs. |
| |
| copyright (c) 1990 Xiaoqiu Huang and Webb Miller |
| |
| For the description of the algorithm, please see the paper |
| "A Time-Efficient, Linear-Space Local Similarity Algorithm" |
| (to appear in Advances in Applied Mathematics) |
| Xiaoqiu Huang and Webb Miller |
| Department of Computer Science |
| The Pennsylvania State University |
| University Park, PA 16802 |
| |
| An early version of this program is presented in the paper |
| "A Space-Efficient Algorithm for Local Similarities" |
| (to appear in Computer Applications in the Biosciences) |
| Xiaoqiu Huang, Ross C. Hardison and Webb Miller |
| Department of Computer Science |
| Department of Molecular and Cell Biology |
| The Pennsylvania State University |
| University Park, PA 16802 |
| |
| The SIM program finds k best non-intersecting alignments between |
| two sequences or within one sequence. Using dynamic programming |
| techniques, SIM is guaranteed to find optimal alignments. The |
| alignments are reported in order of similarity score, with the |
| highest scoring alignment first. The k best alignments share no |
| aligned pairs. SIM requires space proportional to the sum of the |
| input sequence lengths and the output alignment lengths. Thus |
| SIM can handle sequences of tens of thousands, even a hundred of |
| thousands, of base pairs on a workstation. |
| |
| Users can supply values for the parameters: |
| M = cost of a matching aligned pair (default = 1) |
| I = cost of a transition (default is -1) |
| V = cost of a transversion (default is -1) |
| O = gap open penalty (default is 6.0) |
| E = gap extension penalty (default is 0.2) |
| Thus the score for an N-symbol indel is -(Q + R*N). Values are |
| specified with a command argument like O=5.5, where the given number |
| is rounded by SIM to the nearest tenth. |
| |
| For example, to find 7 best non-intersecting alignments of segments |
| from two sequences in files A and B, and using the above default |
| values except that transistions are scored 0, use the command |
| |
| sim 7 A B I=0 |
| |
| Acknowledgments |
| The functions diff() and display() are from Gene Myers. We made |
| the following modifications: similarity weights (integer), instead of |
| distance weights (float), are used, the aligned pairs already output |
| are not used in the subsequent computation, and the printed alignments |
| have been re-formatted in a number of small ways. |
| */ |
| |
| /***************************************************************/ |
| /* Timer options. You MUST uncomment one of the options below */ |
| /* or compile, for example, with the '-DUNIX' option. */ |
| /***************************************************************/ |
| /* #define Amiga */ |
| /* #define UNIX */ |
| /* #define UNIX_Old */ |
| /* #define VMS */ |
| /* #define BORLAND_C */ |
| /* #define MSC */ |
| /* #define MAC */ |
| /* #define IPSC */ |
| /* #define FORTRAN_SEC */ |
| /* #define GTODay */ |
| /* #define CTimer */ |
| /* #define UXPM */ |
| /* #define MAC_TMgr */ |
| /* #define PARIX */ |
| /* #define POSIX */ |
| /* #define WIN32 */ |
| /***********************/ |
| |
| #include <stdio.h> |
| #include <stdlib.h> |
| #include <string.h> |
| #define long int |
| |
| #ifdef __MINGW32__ |
| #undef UNIX |
| #endif |
| |
| #ifdef __XS1B__ |
| #include <unistd.h> |
| #define CTimer |
| #undef UNIX |
| #endif |
| |
| int SIM(); |
| int big_pass(); |
| int fatal(); |
| int fatalf(); |
| int locate(); |
| int no_cross(); |
| int small_pass(); |
| |
| static const char *mybasename(const char *str) { |
| const char *base = strrchr(str, '/'); |
| return base ? base+1 : str; |
| } |
| |
| /* ckalloc - allocate space; check for success */ |
| #if 1 |
| #define ckalloc(x) malloc(x) |
| #else |
| static char *ckalloc(long amount) { |
| char *malloc(), *p; |
| |
| if ((p = malloc( (unsigned) amount)) == NULL) |
| fatal("Ran out of memory."); |
| return(p); |
| } |
| #endif |
| |
| static char *name1, *name2; /* names of sequence files */ |
| |
| #define CMAX 256 /* for EBCDIC */ |
| /*#define round(x) nint(x) /* round to larger magnitude */ |
| #define round(x) ((x)>0.0 ? (x)+0.5 : (x)-0.5) |
| |
| main(argc, argv) int argc; char *argv[]; |
| { long M, N, K; /* Sequence lengths and k */ |
| char *A, *B; /* Storing two sequences */ |
| int symbol; |
| long V[CMAX][CMAX], Q, R; /* Converted integer weights */ |
| double starttime, benchtime, dtime(); |
| double parm_M, parm_I, parm_V, parm_O, parm_E, v; |
| double atof(); |
| FILE *Bp, *Ap, *Cp, *ckopen(); |
| char *arg1, *arg2, *arg3; |
| |
| if ((Cp = fopen("Output/sim.res","a+")) == NULL) |
| { |
| printf("Can not open Output/sim.res\n\n"); |
| exit(1); |
| } |
| |
| |
| starttime = dtime(); |
| |
| if ( argc < 4) { |
| /*fatal("SIM k file1 file2 [M=] [I=] [V=] [O=] [E=]");*/ |
| arg1 = "8"; |
| arg2 = "tob.38-44"; |
| arg3 = "liv.42-48"; |
| } else { |
| arg1 = argv[1]; |
| arg2 = argv[2]; |
| arg3 = argv[3]; |
| } |
| |
| /* read k: the number of local alignments to find */ |
| sscanf(arg1,"%d", &K); |
| if (K == 0) |
| fatal("specified 0 alignments"); |
| |
| /* determine the sequence lengths */ |
| Ap = ckopen(arg2, "r"); |
| for (M = 0; ( symbol = getc(Ap)) != EOF ; ) |
| if ( symbol != '\n' ) |
| ++M; |
| fclose(Ap); |
| name1 = arg2; |
| |
| /* allocate space for A */ |
| A = ( char * ) ckalloc( (M + 1) * sizeof(char)); |
| |
| /* read the first sequence into A */ |
| Ap = ckopen(arg2, "r"); |
| for (M = 0; ( symbol = getc(Ap)) != EOF ; ) |
| if ( symbol != '\n' ) |
| A[++M] = symbol; |
| |
| if (strcmp(arg2,arg3)) { /* sequences are different */ |
| /* read the second sequence into B */ |
| Bp = ckopen(arg3, "r"); |
| for (N = 0; ( symbol = getc(Bp)) != EOF ; ) |
| if ( symbol != '\n' ) |
| ++N; |
| fclose(Bp); |
| name2 = arg3; |
| B = ( char * ) ckalloc( (N + 1) * sizeof(char)); |
| Bp = ckopen(arg3, "r"); |
| for (N = 0; ( symbol = getc(Bp)) != EOF ; ) |
| if ( symbol != '\n' ) |
| B[++N] = symbol; |
| } |
| |
| parm_M = 1.0; |
| parm_I = -1.0; |
| parm_V = -1.0; |
| parm_O = 6.0; |
| parm_E = 0.2; |
| while (--argc > 3) { |
| if (argv[argc][1] != '=') |
| fatalf("argument %d has improper form", argc); |
| v = atof(argv[argc]+2); |
| switch (argv[argc][0]) { |
| case 'M': parm_M = v; break; |
| case 'I': parm_I = v; break; |
| case 'V': parm_V = v; break; |
| case 'O': parm_O = v; break; |
| case 'E': parm_E = v; break; |
| default: fatal("options are M, I, V, O and E."); |
| } |
| } |
| |
| printf("\t\tSIM output with parameters:\n"); |
| printf("\t\tM = %g, I = %g, V = %g\n\t\tO = %g, E = %g\n\n", |
| parm_M, parm_I, parm_V, parm_O, parm_E); |
| |
| /* set up scoring matrix */ |
| V['A']['A'] = V['C']['C'] = |
| V['G']['G'] = V['T']['T'] = round(10*parm_M); |
| |
| V['A']['G'] = V['G']['A'] = |
| V['C']['T'] = V['T']['C'] = round(10*parm_I); |
| |
| V['A']['C'] = V['A']['T'] = |
| V['C']['A'] = V['C']['G'] = |
| V['G']['C'] = V['G']['T'] = |
| V['T']['A'] = V['T']['G'] = round(10*parm_V); |
| |
| Q = round (10*parm_O); |
| R = round (10*parm_E); |
| |
| if (strcmp(arg2, arg3)) |
| SIM(A,B,M,N,K,V,Q,R,2L); |
| else |
| SIM(A,A,M,M,K,V,Q,R,1L); |
| |
| benchtime = dtime() - starttime; |
| /*printf("\nRun Time (sec) = %9.1lf\n",benchtime);*/ |
| |
| fprintf(Cp,"\n"); |
| fprintf(Cp," Run Line: %s %s %s\n", arg1, arg2, arg3); |
| /*fprintf(Cp," Run Time: %9.1lf (sec)\n",benchtime);*/ |
| fprintf(Cp,"#######################################################\n"); |
| fclose(Cp); |
| exit(0); |
| } |
| |
| static long (*v)[CMAX]; /* substitution scores */ |
| static long q, r; /* gap penalties */ |
| static long qr; /* qr = q + r */ |
| |
| typedef struct ONE { long COL ; struct ONE *NEXT ;} pair, *pairptr; |
| pairptr *row, z; /* for saving used aligned pairs */ |
| static short tt; |
| |
| typedef struct NODE |
| { long SCORE; |
| long STARI; |
| long STARJ; |
| long ENDI; |
| long ENDJ; |
| long TOP; |
| long BOT; |
| long LEFT; |
| long RIGHT; } vertex, *vertexptr; |
| |
| vertexptr *LIST; /* an array for saving k best scores */ |
| vertexptr low = 0; /* lowest score node in LIST */ |
| vertexptr most = 0; /* latestly accessed node in LIST */ |
| static long numnode; /* the number of nodes in LIST */ |
| |
| static long *CC, *DD; /* saving matrix scores */ |
| static long *RR, *SS, *EE, *FF; /* saving start-points */ |
| static long *HH, *WW; /* saving matrix scores */ |
| static long *II, *JJ, *XX, *YY; /* saving start-points */ |
| static long m1, mm, n1, nn; /* boundaries of recomputed area */ |
| static long rl, cl; /* left and top boundaries */ |
| static long min; /* minimum score in LIST */ |
| static short flag; /* indicate if recomputation necessary*/ |
| |
| /* DIAG() assigns value to x if (ii,jj) is never used before */ |
| #define DIAG(ii, jj, x, value) \ |
| { for ( tt = 1, z = row[(ii)]; z != 0; z = z->NEXT ) \ |
| if ( z->COL == (jj) ) \ |
| { tt = 0; break; } \ |
| if ( tt ) \ |
| x = ( value ); \ |
| } |
| |
| /* replace (ss1, xx1, yy1) by (ss2, xx2, yy2) if the latter is large */ |
| #define ORDER(ss1, xx1, yy1, ss2, xx2, yy2) \ |
| { if ( ss1 < ss2 ) \ |
| { ss1 = ss2; xx1 = xx2; yy1 = yy2; } \ |
| else \ |
| if ( ss1 == ss2 ) \ |
| { if ( xx1 < xx2 ) \ |
| { xx1 = xx2; yy1 = yy2; } \ |
| else \ |
| if ( xx1 == xx2 && yy1 < yy2 ) \ |
| yy1 = yy2; \ |
| } \ |
| } |
| |
| /* The following definitions are for function diff() */ |
| |
| long diff(), display(); |
| static long zero = 0; /* long type zero */ |
| |
| #define gap(k) ((k) <= 0 ? 0 : q+r*(k)) /* k-symbol indel score */ |
| |
| static long *sapp; /* Current script append ptr */ |
| static long last; /* Last script op appended */ |
| |
| static long I, J; /* current positions of A ,B */ |
| static long no_mat; /* number of matches */ |
| static long no_mis; /* number of mismatches */ |
| static long al_len; /* length of alignment */ |
| /* Append "Delete k" op */ |
| #define DEL(k) \ |
| { I += k; \ |
| al_len += k; \ |
| if (last < 0) \ |
| last = sapp[-1] -= (k); \ |
| else \ |
| last = *sapp++ = -(k); \ |
| } |
| /* Append "Insert k" op */ |
| #define INS(k) \ |
| { J += k; \ |
| al_len += k; \ |
| if (last < 0) \ |
| { sapp[-1] = (k); *sapp++ = last; } \ |
| else \ |
| last = *sapp++ = (k); \ |
| } |
| |
| /* Append "Replace" op */ |
| #define REP \ |
| { last = *sapp++ = 0; \ |
| al_len += 1; \ |
| } |
| |
| /* SIM(A,B,M,N,K,V,Q,R) reports K best non-intersecting alignments of |
| the segments of A and B in order of similarity scores, where |
| V[a][b] is the score of aligning a and b, and -(Q+R*i) is the score |
| of an i-symbol indel. */ |
| |
| SIM(A,B,M,N,K,V,Q,R,nseq) |
| char A[],B[]; long M,N,K,V[][CMAX],Q,R,nseq; |
| { long endi, endj, stari, starj; /* endpoint and startpoint */ |
| long score; /* the max score in LIST */ |
| long count; /* maximum size of list */ |
| register long i, j; /* row and column indices */ |
| long *S; /* saving operations for diff */ |
| vertexptr cur; /* temporary pointer */ |
| vertexptr findmax(); /* return the largest score node */ |
| |
| /* allocate space for all vectors */ |
| j = (N + 1) * sizeof(long); |
| CC = ( long * ) ckalloc(j); |
| DD = ( long * ) ckalloc(j); |
| RR = ( long * ) ckalloc(j); |
| SS = ( long * ) ckalloc(j); |
| EE = ( long * ) ckalloc(j); |
| FF = ( long * ) ckalloc(j); |
| i = (M + 1) * sizeof(long); |
| HH = ( long * ) ckalloc(i); |
| WW = ( long * ) ckalloc(i); |
| II = ( long * ) ckalloc(i); |
| JJ = ( long * ) ckalloc(i); |
| XX = ( long * ) ckalloc(i); |
| YY = ( long * ) ckalloc(i); |
| S = ( long * ) ckalloc(i + j); |
| row = ( pairptr * ) ckalloc( (M + 1) * sizeof(pairptr)); |
| |
| /* set up list for each row */ |
| for ( i = 1; i <= M; i++ ) |
| if ( nseq == 2 ) |
| row[i] = 0; |
| else |
| { row[i] = z = ( pairptr ) ckalloc( (long) sizeof(pair)); |
| z->COL = i; |
| z->NEXT = 0; |
| } |
| |
| v = V; |
| q = Q; |
| r = R; |
| qr = q + r; |
| |
| LIST = ( vertexptr * ) ckalloc( K * sizeof(vertexptr)); |
| for ( i = 0; i < K ; i++ ) |
| LIST[i] = ( vertexptr ) ckalloc( (long) sizeof(vertex)); |
| |
| printf(" Upper Sequence : %s\n", mybasename(name1)); |
| printf(" Length : %d\n", M); |
| printf(" Lower Sequence : %s\n", mybasename(name2)); |
| printf(" Length : %d\n", N); |
| |
| numnode = min = 0; |
| big_pass(A,B,M,N,K,nseq); |
| |
| /* Report the K best alignments one by one. After each |
| alignment is output, recompute part of the matrix. First |
| determine the size of the area to be recomputed, then do the |
| recomputation */ |
| |
| for ( count = K - 1; count >= 0 ; count-- ) |
| { if ( numnode == 0 ) |
| fatal("The number of alignments computed is too large"); |
| cur = findmax();/* Return a pointer to a node with max score*/ |
| score = cur->SCORE; |
| stari = ++cur->STARI; |
| starj = ++cur->STARJ; |
| endi = cur->ENDI; |
| endj = cur->ENDJ; |
| m1 = cur->TOP; |
| mm = cur->BOT; |
| n1 = cur->LEFT; |
| nn = cur->RIGHT; |
| rl = endi - stari + 1; |
| cl = endj - starj + 1; |
| I = stari - 1; |
| J = starj - 1; |
| sapp = S; |
| last = 0; |
| al_len = 0; |
| no_mat = 0; |
| no_mis = 0; |
| diff(&A[stari]-1, &B[starj]-1,rl,cl,q,q); |
| /* Output the best alignment */ |
| printf("\n*********************************************************\n"); |
| printf(" Number %d Local Alignment\n", K - count); |
| printf(" Similarity Score : %g\n",score/10.0); |
| printf(" Match Percentage : %d%%\n", (100*no_mat)/al_len); |
| printf(" Number of Matches : %d\n", no_mat); |
| printf(" Number of Mismatches : %d\n", no_mis); |
| printf(" Total Length of Gaps : %d\n", al_len-no_mat-no_mis); |
| printf(" Begins at (%d, %d) and Ends at (%d, %d)\n", |
| stari,starj, endi,endj); |
| display(&A[stari]-1,&B[starj]-1,rl,cl,S,stari,starj); |
| fflush(stdout); |
| |
| if ( count ) |
| { flag = 0; |
| locate(A,B,nseq); |
| if ( flag ) |
| small_pass(A,B,count,nseq); |
| } |
| } |
| } |
| |
| /* A big pass to compute K best classes */ |
| |
| big_pass(A,B,M,N,K,nseq) char A[],B[]; long M,N,K,nseq; |
| { register long i, j; /* row and column indices */ |
| register long c; /* best score at current point */ |
| register long f; /* best score ending with insertion */ |
| register long d; /* best score ending with deletion */ |
| register long p; /* best score at (i-1, j-1) */ |
| register long ci, cj; /* end-point associated with c */ |
| register long di, dj; /* end-point associated with d */ |
| register long fi, fj; /* end-point associated with f */ |
| register long pi, pj; /* end-point associated with p */ |
| long *va; /* pointer to v(A[i], B[j]) */ |
| long addnode(); /* function for inserting a node */ |
| |
| |
| /* Compute the matrix and save the top K best scores in LIST |
| CC : the scores of the current row |
| RR and EE : the starting point that leads to score CC |
| DD : the scores of the current row, ending with deletion |
| SS and FF : the starting point that leads to score DD */ |
| /* Initialize the 0 th row */ |
| for ( j = 1; j <= N ; j++ ) |
| { CC[j] = 0; |
| RR[j] = 0; |
| EE[j] = j; |
| DD[j] = - (q); |
| SS[j] = 0; |
| FF[j] = j; |
| } |
| for ( i = 1; i <= M; i++) |
| { c = 0; /* Initialize column 0 */ |
| f = - (q); |
| ci = fi = i; |
| va = v[A[i]]; |
| if ( nseq == 2 ) |
| { p = 0; |
| pi = i - 1; |
| cj = fj = pj = 0; |
| } |
| else |
| { p = CC[i]; |
| pi = RR[i]; |
| pj = EE[i]; |
| cj = fj = i; |
| } |
| for ( j = (nseq == 2 ? 1 : (i+1)) ; j <= N ; j++ ) |
| { f = f - r; |
| c = c - qr; |
| ORDER(f, fi, fj, c, ci, cj) |
| c = CC[j] - qr; |
| ci = RR[j]; |
| cj = EE[j]; |
| d = DD[j] - r; |
| di = SS[j]; |
| dj = FF[j]; |
| ORDER(d, di, dj, c, ci, cj) |
| c = 0; |
| DIAG(i, j, c, p+va[B[j]]) /* diagonal */ |
| if ( c <= 0 ) |
| { c = 0; ci = i; cj = j; } |
| else |
| { ci = pi; cj = pj; } |
| ORDER(c, ci, cj, d, di, dj) |
| ORDER(c, ci, cj, f, fi, fj) |
| p = CC[j]; |
| CC[j] = c; |
| pi = RR[j]; |
| pj = EE[j]; |
| RR[j] = ci; |
| EE[j] = cj; |
| DD[j] = d; |
| SS[j] = di; |
| FF[j] = dj; |
| if ( c > min ) /* add the score into list */ |
| min = addnode(c, ci, cj, i, j, K, min); |
| } |
| } |
| } |
| |
| /* Determine the left and top boundaries of the recomputed area */ |
| |
| locate(A,B,nseq) char A[],B[]; long nseq; |
| { register long i, j; /* row and column indices */ |
| register long c; /* best score at current point */ |
| register long f; /* best score ending with insertion */ |
| register long d; /* best score ending with deletion */ |
| register long p; /* best score at (i-1, j-1) */ |
| register long ci, cj; /* end-point associated with c */ |
| register long di, dj; /* end-point associated with d */ |
| register long fi, fj; /* end-point associated with f */ |
| register long pi, pj; /* end-point associated with p */ |
| short cflag, rflag; /* for recomputation */ |
| long *va; /* pointer to v(A[i], B[j]) */ |
| long addnode(); /* function for inserting a node */ |
| long limit; /* the bound on j */ |
| |
| /* Reverse pass |
| rows |
| CC : the scores on the current row |
| RR and EE : the endpoints that lead to CC |
| DD : the deletion scores |
| SS and FF : the endpoints that lead to DD |
| |
| columns |
| HH : the scores on the current columns |
| II and JJ : the endpoints that lead to HH |
| WW : the deletion scores |
| XX and YY : the endpoints that lead to WW |
| */ |
| for ( j = nn; j >= n1 ; j-- ) |
| { CC[j] = 0; |
| EE[j] = j; |
| DD[j] = - (q); |
| FF[j] = j; |
| if ( nseq == 2 || j > mm ) |
| RR[j] = SS[j] = mm + 1; |
| else |
| RR[j] = SS[j] = j; |
| } |
| |
| for ( i = mm; i >= m1; i-- ) |
| { c = p = 0; |
| f = - (q); |
| ci = fi = i; |
| pi = i + 1; |
| cj = fj = pj = nn + 1; |
| va = v[A[i]]; |
| if ( nseq == 2 || n1 > i ) |
| limit = n1; |
| else |
| limit = i + 1; |
| for ( j = nn; j >= limit ; j-- ) |
| { f = f - r; |
| c = c - qr; |
| ORDER(f, fi, fj, c, ci, cj) |
| c = CC[j] - qr; |
| ci = RR[j]; |
| cj = EE[j]; |
| d = DD[j] - r; |
| di = SS[j]; |
| dj = FF[j]; |
| ORDER(d, di, dj, c, ci, cj) |
| c = 0; |
| DIAG(i, j, c, p+va[B[j]]) /* diagonal */ |
| if ( c <= 0 ) |
| { c = 0; ci = i; cj = j; } |
| else |
| { ci = pi; cj = pj; } |
| ORDER(c, ci, cj, d, di, dj) |
| ORDER(c, ci, cj, f, fi, fj) |
| p = CC[j]; |
| CC[j] = c; |
| pi = RR[j]; |
| pj = EE[j]; |
| RR[j] = ci; |
| EE[j] = cj; |
| DD[j] = d; |
| SS[j] = di; |
| FF[j] = dj; |
| if ( c > min ) |
| flag = 1; |
| } |
| if ( nseq == 2 || i < n1 ) |
| { HH[i] = CC[n1]; |
| II[i] = RR[n1]; |
| JJ[i] = EE[n1]; |
| WW[i] = DD[n1]; |
| XX[i] = SS[n1]; |
| YY[i] = FF[n1]; |
| } |
| } |
| |
| for ( rl = m1, cl = n1; ; ) |
| { for ( rflag = cflag = 1; (rflag && m1 > 1) || ( cflag && n1 > 1);) |
| { if ( rflag && m1 > 1 ) /* Compute one row */ |
| { rflag = 0; |
| m1--; |
| c = p = 0; |
| f = - (q); |
| ci = fi = m1; |
| pi = m1 + 1; |
| cj = fj = pj = nn + 1; |
| va = v[A[m1]]; |
| for ( j = nn; j >= n1 ; j-- ) |
| { f = f - r; |
| c = c - qr; |
| ORDER(f, fi, fj, c, ci, cj) |
| c = CC[j] - qr; |
| ci = RR[j]; |
| cj = EE[j]; |
| d = DD[j] - r; |
| di = SS[j]; |
| dj = FF[j]; |
| ORDER(d, di, dj, c, ci, cj) |
| c = 0; |
| DIAG(m1, j, c, p+va[B[j]]) /* diagonal */ |
| if ( c <= 0 ) |
| { c = 0; ci = m1; cj = j; } |
| else |
| { ci = pi; cj = pj; } |
| ORDER(c, ci, cj, d, di, dj) |
| ORDER(c, ci, cj, f, fi, fj) |
| p = CC[j]; |
| CC[j] = c; |
| pi = RR[j]; |
| pj = EE[j]; |
| RR[j] = ci; |
| EE[j] = cj; |
| DD[j] = d; |
| SS[j] = di; |
| FF[j] = dj; |
| if ( c > min ) |
| flag = 1; |
| if ( ! rflag && (ci > rl && cj > cl || di > rl && dj>cl |
| || fi > rl && fj > cl)) |
| rflag = 1; |
| } |
| HH[m1] = CC[n1]; |
| II[m1] = RR[n1]; |
| JJ[m1] = EE[n1]; |
| WW[m1] = DD[n1]; |
| XX[m1] = SS[n1]; |
| YY[m1] = FF[n1]; |
| if ( ! cflag && ( ci > rl && cj > cl || di > rl && dj > cl |
| || fi > rl && fj > cl ) ) |
| cflag = 1; |
| } |
| |
| if ( nseq == 1 && n1 == (m1 + 1) && ! rflag ) |
| cflag = 0; |
| if ( cflag && n1 > 1 ) /* Compute one column */ |
| { cflag = 0; |
| n1--; |
| c = 0; |
| f = - (q); |
| cj = fj = n1; |
| va = v[B[n1]]; |
| if ( nseq == 2 || mm < n1 ) |
| { p = 0; |
| ci = fi = pi = mm + 1; |
| pj = n1 + 1; |
| limit = mm; |
| } |
| else |
| { p = HH[n1]; |
| pi = II[n1]; |
| pj = JJ[n1]; |
| ci = fi = n1; |
| limit = n1 - 1; |
| } |
| for ( i = limit; i >= m1 ; i-- ) |
| { f = f - r; |
| c = c - qr; |
| ORDER(f, fi, fj, c, ci, cj) |
| c = HH[i] - qr; |
| ci = II[i]; |
| cj = JJ[i]; |
| d = WW[i] - r; |
| di = XX[i]; |
| dj = YY[i]; |
| ORDER(d, di, dj, c, ci, cj) |
| c = 0; |
| DIAG(i, n1, c, p+va[A[i]]) |
| if ( c <= 0 ) |
| { c = 0; ci = i; cj = n1; } |
| else |
| { ci = pi; cj = pj; } |
| ORDER(c, ci, cj, d, di, dj) |
| ORDER(c, ci, cj, f, fi, fj) |
| p = HH[i]; |
| HH[i] = c; |
| pi = II[i]; |
| pj = JJ[i]; |
| II[i] = ci; |
| JJ[i] = cj; |
| WW[i] = d; |
| XX[i] = di; |
| YY[i] = dj; |
| if ( c > min ) |
| flag = 1; |
| if (! cflag && (ci > rl && cj > cl || di > rl && dj >cl |
| || fi > rl && fj > cl ) ) |
| cflag = 1; |
| } |
| CC[n1] = HH[m1]; |
| RR[n1] = II[m1]; |
| EE[n1] = JJ[m1]; |
| DD[n1] = WW[m1]; |
| SS[n1] = XX[m1]; |
| FF[n1] = YY[m1]; |
| if ( ! rflag && ( ci > rl && cj > cl || di > rl && dj > cl |
| || fi > rl && fj > cl )) |
| rflag = 1; |
| } |
| } |
| if ( m1 == 1 && n1 == 1 || no_cross() ) |
| break; |
| } |
| m1--; |
| n1--; |
| } |
| |
| /* recompute the area on forward pass */ |
| small_pass(A,B,count,nseq) char A[], B[]; long count, nseq; |
| { register long i, j; /* row and column indices */ |
| register long c; /* best score at current point */ |
| register long f; /* best score ending with insertion */ |
| register long d; /* best score ending with deletion */ |
| register long p; /* best score at (i-1, j-1) */ |
| register long ci, cj; /* end-point associated with c */ |
| register long di, dj; /* end-point associated with d */ |
| register long fi, fj; /* end-point associated with f */ |
| register long pi, pj; /* end-point associated with p */ |
| long *va; /* pointer to v(A[i], B[j]) */ |
| long addnode(); /* function for inserting a node */ |
| long limit; /* lower bound on j */ |
| |
| for ( j = n1 + 1; j <= nn ; j++ ) |
| { CC[j] = 0; |
| RR[j] = m1; |
| EE[j] = j; |
| DD[j] = - (q); |
| SS[j] = m1; |
| FF[j] = j; |
| } |
| for ( i = m1 + 1; i <= mm; i++) |
| { c = 0; /* Initialize column 0 */ |
| f = - (q); |
| ci = fi = i; |
| va = v[A[i]]; |
| if ( nseq == 2 || i <= n1 ) |
| { p = 0; |
| pi = i - 1; |
| cj = fj = pj = n1; |
| limit = n1 + 1; |
| } |
| else |
| { p = CC[i]; |
| pi = RR[i]; |
| pj = EE[i]; |
| cj = fj = i; |
| limit = i + 1; |
| } |
| for ( j = limit ; j <= nn ; j++ ) |
| { f = f - r; |
| c = c - qr; |
| ORDER(f, fi, fj, c, ci, cj) |
| c = CC[j] - qr; |
| ci = RR[j]; |
| cj = EE[j]; |
| d = DD[j] - r; |
| di = SS[j]; |
| dj = FF[j]; |
| ORDER(d, di, dj, c, ci, cj) |
| c = 0; |
| DIAG(i, j, c, p+va[B[j]]) /* diagonal */ |
| if ( c <= 0 ) |
| { c = 0; ci = i; cj = j; } |
| else |
| { ci = pi; cj = pj; } |
| ORDER(c, ci, cj, d, di, dj) |
| ORDER(c, ci, cj, f, fi, fj) |
| p = CC[j]; |
| CC[j] = c; |
| pi = RR[j]; |
| pj = EE[j]; |
| RR[j] = ci; |
| EE[j] = cj; |
| DD[j] = d; |
| SS[j] = di; |
| FF[j] = dj; |
| if ( c > min ) /* add the score into list */ |
| min = addnode(c, ci, cj, i, j, count, min); |
| } |
| } |
| } |
| |
| /* Add a new node into list. */ |
| |
| long addnode(c, ci, cj, i, j, K, cost) long c, ci, cj, i, j, K, cost; |
| { short found; /* 1 if the node is in LIST */ |
| register long d; |
| |
| found = 0; |
| if ( most != 0 && most->STARI == ci && most->STARJ == cj ) |
| found = 1; |
| else |
| for ( d = 0; d < numnode ; d++ ) |
| { most = LIST[d]; |
| if ( most->STARI == ci && most->STARJ == cj ) |
| { found = 1; |
| break; |
| } |
| } |
| if ( found ) |
| { if ( most->SCORE < c ) |
| { most->SCORE = c; |
| most->ENDI = i; |
| most->ENDJ = j; |
| } |
| if ( most->TOP > i ) most->TOP = i; |
| if ( most->BOT < i ) most->BOT = i; |
| if ( most->LEFT > j ) most->LEFT = j; |
| if ( most->RIGHT < j ) most->RIGHT = j; |
| } |
| else |
| { if ( numnode == K ) /* list full */ |
| most = low; |
| else |
| most = LIST[numnode++]; |
| most->SCORE = c; |
| most->STARI = ci; |
| most->STARJ = cj; |
| most->ENDI = i; |
| most->ENDJ = j; |
| most->TOP = most->BOT = i; |
| most->LEFT = most->RIGHT = j; |
| } |
| if ( numnode == K ) |
| { if ( low == most || ! low ) |
| { for ( low = LIST[0], d = 1; d < numnode ; d++ ) |
| if ( LIST[d]->SCORE < low->SCORE ) |
| low = LIST[d]; |
| } |
| return ( low->SCORE ) ; |
| } |
| else |
| return cost; |
| } |
| |
| /* Find and remove the largest score in list */ |
| |
| vertexptr findmax() |
| { vertexptr cur; |
| register long i, j; |
| |
| for ( j = 0, i = 1; i < numnode ; i++ ) |
| if ( LIST[i]->SCORE > LIST[j]->SCORE ) |
| j = i; |
| cur = LIST[j]; |
| if ( j != --numnode ) |
| { LIST[j] = LIST[numnode]; |
| LIST[numnode] = cur; |
| } |
| most = LIST[0]; |
| if ( low == cur ) low = LIST[0]; |
| return ( cur ); |
| } |
| |
| /* return 1 if no node in LIST share vertices with the area */ |
| |
| no_cross() |
| { vertexptr cur; |
| register long i; |
| |
| for ( i = 0; i < numnode; i++ ) |
| { cur = LIST[i]; |
| if ( cur->STARI <= mm && cur->STARJ <= nn && cur->BOT >= m1-1 && |
| cur->RIGHT >= n1-1 && (cur->STARI < rl || cur->STARJ < cl)) |
| { if ( cur->STARI < rl ) rl = cur->STARI; |
| if ( cur->STARJ < cl ) cl = cur->STARJ; |
| flag = 1; |
| break; |
| } |
| } |
| if ( i == numnode ) |
| return 1; |
| else |
| return 0; |
| } |
| |
| /* diff(A,B,M,N,tb,te) returns the score of an optimum conversion between |
| A[1..M] and B[1..N] that begins(ends) with a delete if tb(te) is zero |
| and appends such a conversion to the current script. */ |
| |
| long diff(A,B,M,N,tb,te) char *A, *B; long M, N; long tb, te; |
| |
| { long midi, midj, type; /* Midpoint, type, and cost */ |
| long midc; |
| |
| { register long i, j; |
| register long c, e, d, s; |
| long t, *va; |
| |
| /* Boundary cases: M <= 1 or N == 0 */ |
| |
| if (N <= 0) |
| { if (M > 0) DEL(M) |
| return - gap(M); |
| } |
| if (M <= 1) |
| { if (M <= 0) |
| { INS(N); |
| return - gap(N); |
| } |
| if (tb > te) tb = te; |
| midc = - (tb + r + gap(N) ); |
| midj = 0; |
| va = v[A[1]]; |
| for (j = 1; j <= N; j++) |
| { for ( tt = 1, z = row[I+1]; z != 0; z = z->NEXT ) |
| if ( z->COL == j+J ) |
| { tt = 0; break; } |
| if ( tt ) |
| { c = va[B[j]] - ( gap(j-1) + gap(N-j) ); |
| if (c > midc) |
| { midc = c; |
| midj = j; |
| } |
| } |
| } |
| if (midj == 0) |
| { INS(N) DEL(1) } |
| else |
| { if (midj > 1) INS(midj-1) |
| REP |
| if ( A[1] == B[midj] ) |
| no_mat += 1; |
| else |
| no_mis += 1; |
| /* mark (A[I],B[J]) as used: put J into list row[I] */ |
| I++; J++; |
| z = ( pairptr ) ckalloc( (long) sizeof(pair)); |
| z->COL = J; |
| z->NEXT = row[I]; |
| row[I] = z; |
| if (midj < N) INS(N-midj) |
| } |
| return midc; |
| } |
| |
| /* Divide: Find optimum midpoint (midi,midj) of cost midc */ |
| |
| midi = M/2; /* Forward phase: */ |
| CC[0] = 0; /* Compute C(M/2,k) & D(M/2,k) for all k */ |
| t = -q; |
| for (j = 1; j <= N; j++) |
| { CC[j] = t = t-r; |
| DD[j] = t-q; |
| } |
| t = -tb; |
| for (i = 1; i <= midi; i++) |
| { s = CC[0]; |
| CC[0] = c = t = t-r; |
| e = t-q; |
| va = v[A[i]]; |
| for (j = 1; j <= N; j++) |
| { if ((c = c - qr) > (e = e - r)) e = c; |
| if ((c = CC[j] - qr) > (d = DD[j] - r)) d = c; |
| DIAG(i+I, j+J, c, s+va[B[j]]) |
| if (c < d) c = d; |
| if (c < e) c = e; |
| s = CC[j]; |
| CC[j] = c; |
| DD[j] = d; |
| } |
| } |
| DD[0] = CC[0]; |
| |
| RR[N] = 0; /* Reverse phase: */ |
| t = -q; /* Compute R(M/2,k) & S(M/2,k) for all k */ |
| for (j = N-1; j >= 0; j--) |
| { RR[j] = t = t-r; |
| SS[j] = t-q; |
| } |
| t = -te; |
| for (i = M-1; i >= midi; i--) |
| { s = RR[N]; |
| RR[N] = c = t = t-r; |
| e = t-q; |
| va = v[A[i+1]]; |
| for (j = N-1; j >= 0; j--) |
| { if ((c = c - qr) > (e = e - r)) e = c; |
| if ((c = RR[j] - qr) > (d = SS[j] - r)) d = c; |
| DIAG(i+1+I, j+1+J, c, s+va[B[j+1]]) |
| if (c < d) c = d; |
| if (c < e) c = e; |
| s = RR[j]; |
| RR[j] = c; |
| SS[j] = d; |
| } |
| } |
| SS[N] = RR[N]; |
| |
| midc = CC[0]+RR[0]; /* Find optimal midpoint */ |
| midj = 0; |
| type = 1; |
| for (j = 0; j <= N; j++) |
| if ((c = CC[j] + RR[j]) >= midc) |
| if (c > midc || CC[j] != DD[j] && RR[j] == SS[j]) |
| { midc = c; |
| midj = j; |
| } |
| for (j = N; j >= 0; j--) |
| if ((c = DD[j] + SS[j] + q) > midc) |
| { midc = c; |
| midj = j; |
| type = 2; |
| } |
| } |
| |
| /* Conquer: recursively around midpoint */ |
| |
| if (type == 1) |
| { diff(A,B,midi,midj,tb,q); |
| diff(A+midi,B+midj,M-midi,N-midj,q,te); |
| } |
| else |
| { diff(A,B,midi-1,midj,tb,zero); |
| DEL(2); |
| diff(A+midi+1,B+midj,M-midi-1,N-midj,zero,te); |
| } |
| return midc; |
| } |
| |
| /* Alignment display routine */ |
| |
| static char ALINE[51], BLINE[51], CLINE[51]; |
| |
| long display(A,B,M,N,S,AP,BP) char A[], B[]; long M, N; long S[], AP, BP; |
| { register char *a, *b, *c; |
| register long i, j, op; |
| long lines, ap, bp; |
| char x, y, z; |
| |
| |
| i = j = op = lines = 0; |
| ap = AP; |
| bp = BP; |
| a = ALINE; |
| b = BLINE; |
| c = CLINE; |
| while (i < M || j < N) |
| { if (op == 0 && *S == 0) |
| { op = *S++; |
| *a = A[++i]; |
| *b = B[++j]; |
| x = *a++; |
| y = *b++; |
| z = ' '; |
| if (x == y) |
| z = '|'; |
| if ( |
| (x == 'A' && y == 'G') || |
| (x == 'C' && y == 'T') || |
| (x == 'G' && y == 'A') || |
| (x == 'T' && y == 'C') |
| ) z = ':'; |
| *c++ = z; |
| } |
| else |
| { if (op == 0) |
| op = *S++; |
| if (op > 0) |
| { *a++ = ' '; |
| *b++ = B[++j]; |
| op--; |
| } |
| else |
| { *a++ = A[++i]; |
| *b++ = ' '; |
| op++; |
| } |
| *c++ = '-'; |
| } |
| if (a >= ALINE+50 || i >= M && j >= N) |
| { *a = *b = *c = '\0'; |
| printf("\n%5d ",50*lines++); |
| for (b = ALINE+10; b <= a; b += 10) |
| printf(" . :"); |
| if (b <= a+5) |
| printf(" ."); |
| printf("\n%5d %s\n %s\n%5d %s\n",ap,ALINE,CLINE,bp,BLINE); |
| ap = AP + i; |
| bp = BP + j; |
| a = ALINE; |
| b = BLINE; |
| c = CLINE; |
| } |
| } |
| } |
| |
| /* lib.c - library of C procedures. */ |
| |
| /* fatal - print message and die */ |
| fatal(msg) |
| char *msg; |
| { |
| fprintf(stderr, "%s\n", msg); |
| exit(1); |
| } |
| |
| /* fatalf - format message, print it, and die */ |
| fatalf(msg, val) |
| char *msg, *val; |
| { |
| fprintf(stderr, msg, val); |
| putc('\n', stderr); |
| exit(1); |
| } |
| |
| /* ckopen - open file; check for success */ |
| FILE *ckopen(name, mode) |
| char *name, *mode; |
| { |
| FILE *fopen(), *fp; |
| |
| if ((fp = fopen(name, mode)) == NULL) |
| fatalf("Cannot open %s.", name); |
| return(fp); |
| } |
| |
| |
| /*****************************************************/ |
| /* Various timer routines. */ |
| /* Al Aburto, aburto@nosc.mil, 08 Oct 1996 */ |
| /* */ |
| /* t = dtime() outputs the current time in seconds. */ |
| /* Use CAUTION as some of these routines will mess */ |
| /* up when timing across the hour mark!!! */ |
| /* */ |
| /* For timing I use the 'user' time whenever */ |
| /* possible. Using 'user+sys' time is a separate */ |
| /* issue. */ |
| /* */ |
| /* Example Usage: */ |
| /* [timer options added here] */ |
| /* main() */ |
| /* { */ |
| /* double starttime,benchtime,dtime(); */ |
| /* */ |
| /* starttime = dtime(); */ |
| /* [routine to time] */ |
| /* benchtime = dtime() - starttime; */ |
| /* } */ |
| /* */ |
| /* [timer code below added here] */ |
| /*****************************************************/ |
| |
| /*********************************/ |
| /* Timer code. */ |
| /*********************************/ |
| /*******************/ |
| /* Amiga dtime() */ |
| /*******************/ |
| #ifdef Amiga |
| #include <ctype.h> |
| #define HZ 50 |
| |
| double dtime() |
| { |
| double q; |
| |
| struct tt { |
| long days; |
| long minutes; |
| long ticks; |
| } tt; |
| |
| DateStamp(&tt); |
| |
| q = ((double)(tt.ticks + (tt.minutes * 60L * 50L))) / (double)HZ; |
| |
| return q; |
| } |
| #endif |
| #undef long |
| /*****************************************************/ |
| /* UNIX dtime(). This is the preferred UNIX timer. */ |
| /* Provided by: Markku Kolkka, mk59200@cc.tut.fi */ |
| /* HP-UX Addition by: Bo Thide', bt@irfu.se */ |
| /*****************************************************/ |
| #ifdef UNIX |
| #include <sys/time.h> |
| #include <sys/resource.h> |
| |
| #ifdef hpux |
| #include <sys/syscall.h> |
| #define getrusage(a,b) syscall(SYS_getrusage,a,b) |
| #endif |
| |
| struct rusage rusage; |
| |
| double dtime() |
| { |
| double q; |
| |
| getrusage(RUSAGE_SELF,&rusage); |
| |
| q = (double)(rusage.ru_utime.tv_sec); |
| q = q + (double)(rusage.ru_utime.tv_usec) * 1.0e-06; |
| |
| return q; |
| } |
| #endif |
| |
| /***************************************************/ |
| /* UNIX_Old dtime(). This is the old UNIX timer. */ |
| /* Use only if absolutely necessary as HZ may be */ |
| /* ill defined on your system. */ |
| /***************************************************/ |
| #ifdef UNIX_Old |
| #include <sys/types.h> |
| #include <sys/times.h> |
| #include <sys/param.h> |
| |
| #ifndef HZ |
| #define HZ 60 |
| #endif |
| |
| struct tms tms; |
| |
| double dtime() |
| { |
| double q; |
| |
| times(&tms); |
| |
| q = (double)(tms.tms_utime) / (double)HZ; |
| |
| return q; |
| } |
| #endif |
| |
| /*********************************************************/ |
| /* VMS dtime() for VMS systems. */ |
| /* Provided by: RAMO@uvphys.phys.UVic.CA */ |
| /* Some people have run into problems with this timer. */ |
| /*********************************************************/ |
| #ifdef VMS |
| #include time |
| |
| #ifndef HZ |
| #define HZ 100 |
| #endif |
| |
| struct tbuffer_t |
| { |
| int proc_user_time; |
| int proc_system_time; |
| int child_user_time; |
| int child_system_time; |
| }; |
| struct tbuffer_t tms; |
| |
| double dtime() |
| { |
| double q; |
| |
| times(&tms); |
| |
| q = (double)(tms.proc_user_time) / (double)HZ; |
| |
| return q; |
| } |
| #endif |
| |
| /******************************/ |
| /* BORLAND C dtime() for DOS */ |
| /******************************/ |
| #ifdef BORLAND_C |
| #include <ctype.h> |
| #include <dos.h> |
| #include <time.h> |
| |
| #define HZ 100 |
| struct time tnow; |
| |
| double dtime() |
| { |
| double q; |
| |
| gettime(&tnow); |
| |
| q = 60.0 * (double)(tnow.ti_min); |
| q = q + (double)(tnow.ti_sec); |
| q = q + (double)(tnow.ti_hund)/(double)HZ; |
| |
| return q; |
| } |
| #endif |
| |
| /**************************************/ |
| /* Microsoft C (MSC) dtime() for DOS */ |
| /**************************************/ |
| #ifdef MSC |
| #include <time.h> |
| #include <ctype.h> |
| |
| #define HZ CLOCKS_PER_SEC |
| clock_t tnow; |
| |
| double dtime() |
| { |
| double q; |
| |
| tnow = clock(); |
| |
| q = (double)tnow / (double)HZ; |
| |
| return q; |
| } |
| #endif |
| |
| /*************************************/ |
| /* Macintosh (MAC) Think C dtime() */ |
| /*************************************/ |
| #ifdef MAC |
| #include <time.h> |
| |
| #define HZ 60 |
| |
| double dtime() |
| { |
| double q; |
| |
| q = (double)clock() / (double)HZ; |
| |
| return q; |
| } |
| #endif |
| |
| /************************************************************/ |
| /* iPSC/860 (IPSC) dtime() for i860. */ |
| /* Provided by: Dan Yergeau, yergeau@gloworm.Stanford.EDU */ |
| /************************************************************/ |
| #ifdef IPSC |
| extern double dclock(); |
| |
| double dtime() |
| { |
| double q; |
| |
| q = dclock(); |
| |
| return q; |
| } |
| #endif |
| |
| /**************************************************/ |
| /* FORTRAN dtime() for Cray type systems. */ |
| /* This is the preferred timer for Cray systems. */ |
| /**************************************************/ |
| #ifdef FORTRAN_SEC |
| |
| fortran double second(); |
| |
| double dtime() |
| { |
| double q; |
| |
| second(&q); |
| |
| return q; |
| } |
| #endif |
| |
| /***********************************************************/ |
| /* UNICOS C dtime() for Cray UNICOS systems. Don't use */ |
| /* unless absolutely necessary as returned time includes */ |
| /* 'user+system' time. Provided by: R. Mike Dority, */ |
| /* dority@craysea.cray.com */ |
| /***********************************************************/ |
| #ifdef CTimer |
| #include <time.h> |
| |
| double dtime() |
| { |
| double q; |
| clock_t clock(void); |
| |
| q = (double)clock() / (double)CLOCKS_PER_SEC; |
| |
| return q; |
| } |
| #endif |
| |
| /********************************************/ |
| /* Another UNIX timer using gettimeofday(). */ |
| /* However, getrusage() is preferred. */ |
| /********************************************/ |
| #ifdef GTODay |
| #include <sys/time.h> |
| |
| struct timeval tnow; |
| |
| double dtime() |
| { |
| double q; |
| |
| gettimeofday(&tnow,NULL); |
| q = (double)tnow.tv_sec + (double)tnow.tv_usec * 1.0e-6; |
| |
| return q; |
| } |
| #endif |
| |
| /*****************************************************/ |
| /* Fujitsu UXP/M timer. */ |
| /* Provided by: Mathew Lim, ANUSF, M.Lim@anu.edu.au */ |
| /*****************************************************/ |
| #ifdef UXPM |
| #include <sys/types.h> |
| #include <sys/timesu.h> |
| struct tmsu rusage; |
| |
| double dtime() |
| { |
| double q; |
| |
| timesu(&rusage); |
| |
| q = (double)(rusage.tms_utime) * 1.0e-06; |
| |
| return q; |
| } |
| #endif |
| |
| /**********************************************/ |
| /* Macintosh (MAC_TMgr) Think C dtime() */ |
| /* requires Think C Language Extensions or */ |
| /* #include <MacHeaders> in the prefix */ |
| /* provided by Francis H Schiffer 3rd (fhs) */ |
| /* skipschiffer@genie.geis.com */ |
| /**********************************************/ |
| #ifdef MAC_TMgr |
| #include <Timer.h> |
| #include <stdlib.h> |
| |
| static TMTask mgrTimer; |
| static Boolean mgrInited = false; |
| static double mgrClock; |
| |
| #define RMV_TIMER RmvTime( (QElemPtr)&mgrTimer ) |
| #define MAX_TIME 1800000000L |
| /* MAX_TIME limits time between calls to */ |
| /* dtime( ) to no more than 30 minutes */ |
| /* this limitation could be removed by */ |
| /* creating a completion routine to sum */ |
| /* 30 minute segments (fhs 1994 feb 9) */ |
| |
| static void Remove_timer( ) |
| { |
| RMV_TIMER; |
| mgrInited = false; |
| } |
| double dtime( ) |
| { |
| if( mgrInited ) { |
| RMV_TIMER; |
| mgrClock += (MAX_TIME + mgrTimer.tmCount)*1.0e-6; |
| } else { |
| if( _atexit( &Remove_timer ) == 0 ) mgrInited = true; |
| mgrClock = 0.0; |
| } |
| if( mgrInited ) { |
| mgrTimer.tmAddr = NULL; |
| mgrTimer.tmCount = 0; |
| mgrTimer.tmWakeUp = 0; |
| mgrTimer.tmReserved = 0; |
| InsTime( (QElemPtr)&mgrTimer ); |
| PrimeTime( (QElemPtr)&mgrTimer, -MAX_TIME ); |
| } |
| return( mgrClock ); |
| } |
| #endif |
| |
| /***********************************************************/ |
| /* Parsytec GCel timer. */ |
| /* Provided by: Georg Wambach, gw@informatik.uni-koeln.de */ |
| /***********************************************************/ |
| #ifdef PARIX |
| #include <sys/time.h> |
| |
| double dtime() |
| { |
| double q; |
| |
| q = (double) (TimeNowHigh()) / (double) CLK_TCK_HIGH; |
| |
| return q; |
| } |
| #endif |
| |
| /************************************************/ |
| /* Sun Solaris POSIX dtime() routine */ |
| /* Provided by: Case Larsen, CTLarsen.lbl.gov */ |
| /************************************************/ |
| #ifdef POSIX |
| #include <sys/time.h> |
| #include <sys/resource.h> |
| #include <sys/rusage.h> |
| |
| #ifdef __hpux |
| #include <sys/syscall.h> |
| #endif |
| |
| struct rusage rusage; |
| |
| double dtime() |
| { |
| double q; |
| |
| getrusage(RUSAGE_SELF,&rusage); |
| |
| q = (double)(rusage.ru_utime.tv_sec); |
| q = q + (double)(rusage.ru_utime.tv_nsec) * 1.0e-09; |
| |
| return q; |
| } |
| #endif |
| |
| |
| /****************************************************/ |
| /* Windows NT (32 bit) dtime() routine */ |
| /* Provided by: Piers Haken, piersh@microsoft.com */ |
| /****************************************************/ |
| #ifdef WIN32 |
| #include <windows.h> |
| |
| double dtime(void) |
| { |
| double q; |
| |
| q = (double)GetTickCount() * 1.0e-03; |
| |
| return q; |
| } |
| #endif |
| |
| /*------ Say goodnight Mary! ------*/ |
| |
