1 /* bestscor.c 13-Mar-1985 */
2
3 /* copyright (C) 1983 William R. Pearson */
4
5 #include <stdio.h>
6 #include <stdlib.h>
7 #include <string.h>
8 #include <ctype.h>
9
10 #define TRUE 1
11 #define FALSE 0
12
13 #define MAXTST 1000 /* longest test sequence */
14 #define MAXLIB 3000 /* longest library sequence (not used) */
15 #define MAXDIAG 4000 /* sum of test and library sequence */
16
17 char *aa0;
18 int n0;
19 long sq0off=1;
20
21 #include "upam.gbl"
22 #define XTERNAL
23 #include "uascii.gbl"
24 int histint=2;
25 int bestscale=200;
26 int bkfact=5;
27 int scfact=4;
28 int bktup=2;
29 int ktmax=2;
30 int bestmax=50;
31 int bestoff=27; /* values for calculating bestcut */
32 int dnaseq = 0;
33
34 extern int optind;
35 char smstr[40], *smptr;
36
main(argc,argv)37 main(argc, argv)
38 int argc; char **argv;
39 {
40 char rline[40], tname[40];
41
42 initenv(argc,argv);
43
44 if (argc-optind < 2) {
45 printf(" bestscor calculates the score of a 100%% identical match\n");
46 printf(" test sequence file name: ");
47 fgets(tname,sizeof(tname),stdin);
48 if (tname[strlen(tname)-1]=='\n') tname[strlen(tname)-1]='\0';
49 }
50 else {
51 strncpy(tname,argv[optind+1],sizeof(tname));
52 }
53
54 if ((aa0=calloc((size_t)MAXDIAG,sizeof(char)))==0) {
55 printf(" cannot allocate sequence array\n");
56 exit(1);
57 }
58
59 if ((n0=getseq(tname,aa0,MAXDIAG,&dnaseq))==0) {
60 printf(" %s : %s sequence not found\n",tname,sqtype);
61 exit(1);
62 }
63
64 resetp(dnaseq);
65
66 initpam2(); /* convert 1-d pam to 2-d pam2 */
67
68 printf(" %s : %4d %s\n",tname, n0,sqnam);
69 printf(" 100%% identical score using %s is %d\n",
70 smptr,shscore(aa0,n0));
71 }
72
73 extern int *sascii, nascii[], aascii[];
74
initenv(argc,argv)75 initenv(argc,argv)
76 int argc; char **argv;
77 {
78 char *cptr, *getenv();
79 int copt, getopt();
80 extern char *optarg;
81
82 sascii = aascii;
83 pam = abl50;
84 strncpy(smstr,"BLOSUM50",sizeof(smstr));
85 smptr=smstr;
86 sq = aa;
87 hsq = haa;
88 nsq = naa;
89 dnaseq = 0;
90
91 while ((copt=getopt(argc,argv,"s:"))!=EOF)
92 switch(copt) {
93 case 's': strncpy(smstr,optarg,sizeof(smstr));
94 smptr=smstr;
95 if (initpam(smptr)) dnaseq= -1;
96 else smptr="\0";
97 break;
98 default : fprintf(stderr," illegal option -%c\n",copt);
99 }
100 optind--;
101
102 if (dnaseq>=0) {
103 if ((smptr=getenv("SMATRIX"))!=NULL && initpam(smptr))
104 dnaseq = -1;
105 else {
106 if (dnaseq == 0 ) smptr=smstr;
107 else smptr="DNA";
108 }
109 }
110
111
112 if (strlen(smptr)>0) fprintf(stderr," using matrix file %s\n",smptr);
113 }
114
resetp(dnaseq)115 resetp(dnaseq)
116 int dnaseq;
117 {
118 if (dnaseq==1) {
119 pam = npam;
120 strncpy(smstr,"DNA",sizeof(smstr));
121 smptr = smstr;
122 }
123 }
124
shscore(aa0,n0)125 shscore(aa0,n0) /* calculate the 100% identical score */
126 char *aa0; int n0;
127 {
128 int i, sum;
129 for (i=0,sum=0; i<n0; i++)
130 sum += pam2[aa0[i]][aa0[i]];
131 return sum;
132 }
133
initpam2()134 initpam2()
135 {
136 int i, j, k;
137
138 k=0;
139 for (i=0; i<nsq; i++)
140 for (j=0; j<=i; j++)
141 pam2[j][i] = pam2[i][j] = pam[k++];
142 }
143
144