--- /dev/null
+
+/* Program COILS version 2.1 */
+ /* written by J.Lupas, 22 JUN 1993 */
+ /* edited on 15 JAN 1994 */
+ /* revised and corrected 4 APR 1994 */
+ /* incorporates the option to use either the old MTK chart or
+ the new MTIDK chart AND the choice of weighting positions
+ a & d more (2.5 times) */
+ /* 4 output options:
+ - probabilities for window sizes 14, 21 & 28 in columns,
+ - probabilities for window sizes 14, 21 & 28 in rows,
+ - scores only for a user-input window size, or
+ - the probabilities above a user-input cutoff.
+ Requests input and output files from user. */
+ /* transferred to c++ by Larry Harvie, harvie@owl.WPI.EDU */
+ /* adapted to C by Reinhard Schneider, SCHNEIDER@EMBL-Heidelberg.DE */
+
+ /* Adapted by Rob Russell (SB) to make it less user friendly
+ * (UNIX command line options replace command line-driven commands)
+ * Also a FASTA output option: */
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+
+
+ typedef FILE *text;
+ typedef double chartstore[21][8];
+ typedef char residuestore[10001];
+ typedef double calcstore[4][10001];
+ typedef int heptadstore[4][10001];
+
+ int window;
+ int i;
+ char inputfile[80],outputfile[80];
+ text protein,pickup;
+ chartstore chartvalue; /* probability values */
+ residuestore residue; /* stores amino acid sequence */
+ heptadstore heptnum; /* stores heptad position number */
+ calcstore calcnumb; /* stores best scores */
+ int res_total; /* total no. of residues in protein */
+ char junkline[9000]; /* reads opening protein info. */
+ char nameprot[9000];
+ int nomore; /* in case there's junk after last protein */
+ short windows;
+ double winheight;
+ char option,chartchar,weightchar;
+ double ad_weight;
+ double peakmin;
+ int hept_weight;
+ int verbose;
+ calcstore tempcalc;
+ heptadstore temphept;
+ char fr14[10001],fr21[10001],fr28[10001];
+ int pr14[10001],pr21[10001],pr28[10001];
+
+short mod(short x,short y)
+{
+ double x1=(double) x,y1=(double) y;
+
+ return (short) floor(fmod(x,y));
+}
+void Init_newchart()
+ /* store the probability values for each of the 20
+ amino acids to occur in each of the heptad positions */
+{
+chartvalue[1][1]=2.998; chartvalue[1][2]=0.269; chartvalue[1][3]=0.367;
+chartvalue[1][4]=3.852; chartvalue[1][5]=0.510; chartvalue[1][6]=0.514;
+chartvalue[1][7]=0.562;
+chartvalue[2][1]=2.408; chartvalue[2][2]=0.261; chartvalue[2][3]=0.345;
+chartvalue[2][4]=0.931; chartvalue[2][5]=0.402; chartvalue[2][6]=0.440;
+chartvalue[2][7]=0.289;
+chartvalue[3][1]=1.525; chartvalue[3][2]=0.479; chartvalue[3][3]=0.350;
+chartvalue[3][4]=0.887; chartvalue[3][5]=0.286; chartvalue[3][6]=0.350;
+chartvalue[3][7]=0.362;
+chartvalue[4][1]=2.161; chartvalue[4][2]=0.605; chartvalue[4][3]=0.442;
+chartvalue[4][4]=1.441; chartvalue[4][5]=0.607; chartvalue[4][6]=0.457;
+chartvalue[4][7]=0.570;
+chartvalue[5][1]=0.490; chartvalue[5][2]=0.075; chartvalue[5][3]=0.391;
+chartvalue[5][4]=0.639; chartvalue[5][5]=0.125; chartvalue[5][6]=0.081;
+chartvalue[5][7]=0.038;
+chartvalue[6][1]=1.319; chartvalue[6][2]=0.064; chartvalue[6][3]=0.081;
+chartvalue[6][4]=1.526; chartvalue[6][5]=0.204; chartvalue[6][6]=0.118;
+chartvalue[6][7]=0.096;
+chartvalue[7][1]=0.084; chartvalue[7][2]=0.215; chartvalue[7][3]=0.432;
+chartvalue[7][4]=0.111; chartvalue[7][5]=0.153; chartvalue[7][6]=0.367;
+chartvalue[7][7]=0.125;
+chartvalue[8][1]=1.283; chartvalue[8][2]=1.364; chartvalue[8][3]=1.077;
+chartvalue[8][4]=2.219; chartvalue[8][5]=0.490; chartvalue[8][6]=1.265;
+chartvalue[8][7]=0.903;
+chartvalue[9][1]=1.233; chartvalue[9][2]=2.194; chartvalue[9][3]=1.817;
+chartvalue[9][4]=0.611; chartvalue[9][5]=2.095; chartvalue[9][6]=1.686;
+chartvalue[9][7]=2.027;
+chartvalue[10][1]=1.014; chartvalue[10][2]=1.476; chartvalue[10][3]=1.771;
+chartvalue[10][4]=0.114; chartvalue[10][5]=1.667; chartvalue[10][6]=2.006;
+chartvalue[10][7]=1.844;
+chartvalue[11][1]=0.590; chartvalue[11][2]=0.646; chartvalue[11][3]=0.584;
+chartvalue[11][4]=0.842; chartvalue[11][5]=0.307; chartvalue[11][6]=0.611;
+chartvalue[11][7]=0.396;
+chartvalue[12][1]=0.281; chartvalue[12][2]=3.351; chartvalue[12][3]=2.998;
+chartvalue[12][4]=0.789; chartvalue[12][5]=4.868; chartvalue[12][6]=2.735;
+chartvalue[12][7]=3.812;
+chartvalue[13][1]=0.068; chartvalue[13][2]=2.103; chartvalue[13][3]=1.646;
+chartvalue[13][4]=0.182; chartvalue[13][5]=0.664; chartvalue[13][6]=1.581;
+chartvalue[13][7]=1.401;
+chartvalue[14][1]=0.311; chartvalue[14][2]=2.290; chartvalue[14][3]=2.330;
+chartvalue[14][4]=0.811; chartvalue[14][5]=2.596; chartvalue[14][6]=2.155;
+chartvalue[14][7]=2.585;
+chartvalue[15][1]=1.231; chartvalue[15][2]=1.683; chartvalue[15][3]=2.157;
+chartvalue[15][4]=0.197; chartvalue[15][5]=1.653; chartvalue[15][6]=2.430;
+chartvalue[15][7]=2.065;
+chartvalue[16][1]=0.332; chartvalue[16][2]=0.753; chartvalue[16][3]=0.930;
+chartvalue[16][4]=0.424; chartvalue[16][5]=0.734; chartvalue[16][6]=0.801;
+chartvalue[16][7]=0.518;
+chartvalue[17][1]=0.197; chartvalue[17][2]=0.543; chartvalue[17][3]=0.647;
+chartvalue[17][4]=0.680; chartvalue[17][5]=0.905; chartvalue[17][6]=0.643;
+chartvalue[17][7]=0.808;
+chartvalue[18][1]=0.918; chartvalue[18][2]=0.002; chartvalue[18][3]=0.385;
+chartvalue[18][4]=0.440; chartvalue[18][5]=0.138; chartvalue[18][6]=0.432;
+chartvalue[18][7]=0.079;
+chartvalue[19][1]=0.066; chartvalue[19][2]=0.064; chartvalue[19][3]=0.065;
+chartvalue[19][4]=0.747; chartvalue[19][5]=0.006; chartvalue[19][6]=0.115;
+chartvalue[19][7]=0.014;
+chartvalue[20][1]=0.004; chartvalue[20][2]=0.108; chartvalue[20][3]=0.018;
+chartvalue[20][4]=0.006; chartvalue[20][5]=0.010; chartvalue[20][6]=0.004;
+chartvalue[20][7]=0.007;
+}
+void Init_oldchart()
+{
+ chartvalue[1][1]=3.167; chartvalue[1][2]=0.297;
+chartvalue[1][3]=0.398;
+ chartvalue[1][4]=3.902; chartvalue[1][5]=0.585;
+chartvalue[1][6]=0.501;
+ chartvalue[1][7]=0.483;
+ chartvalue[2][1]=2.597; chartvalue[2][2]=0.098;
+chartvalue[2][3]=0.345;
+ chartvalue[2][4]=0.894; chartvalue[2][5]=0.514;
+chartvalue[2][6]=0.471;
+ chartvalue[2][7]=0.431;
+ chartvalue[3][1]=1.665; chartvalue[3][2]=0.403;
+chartvalue[3][3]=0.386;
+ chartvalue[3][4]=0.949; chartvalue[3][5]=0.211;
+chartvalue[3][6]=0.342;
+ chartvalue[3][7]=0.360;
+ chartvalue[4][1]=2.240; chartvalue[4][2]=0.37;
+chartvalue[4][3]=0.480;
+ chartvalue[4][4]=1.409; chartvalue[4][5]=0.541;
+chartvalue[4][6]=0.772;
+ chartvalue[4][7]=0.663;
+ chartvalue[5][1]=0.531; chartvalue[5][2]=0.076;
+chartvalue[5][3]=0.403;
+ chartvalue[5][4]=0.662; chartvalue[5][5]=0.189;
+chartvalue[5][6]=0.106;
+ chartvalue[5][7]=0.013;
+ chartvalue[6][1]=1.417; chartvalue[6][2]=0.090;
+chartvalue[6][3]=0.122;
+ chartvalue[6][4]=1.659; chartvalue[6][5]=0.19; chartvalue[6][6]=0.13;
+ chartvalue[6][7]=0.1550;
+ chartvalue[7][1]=0.045; chartvalue[7][2]=0.275;
+chartvalue[7][3]=0.578;
+ chartvalue[7][4]=0.216; chartvalue[7][5]=0.211;
+chartvalue[7][6]=0.426;
+ chartvalue[7][7]=0.156;
+ chartvalue[8][1]=1.297; chartvalue[8][2]=1.551;
+chartvalue[8][3]=1.084;
+ chartvalue[8][4]=2.612; chartvalue[8][5]=0.377;
+chartvalue[8][6]=1.248;
+ chartvalue[8][7]=0.877;
+ chartvalue[9][1]=1.375; chartvalue[9][2]=2.639;
+chartvalue[9][3]=1.763;
+ chartvalue[9][4]=0.191; chartvalue[9][5]=1.815;
+chartvalue[9][6]=1.961;
+ chartvalue[9][7]=2.795;
+ chartvalue[10][1]=0.659; chartvalue[10][2]=1.163;
+chartvalue[10][3]=1.210;
+ chartvalue[10][4]=0.031; chartvalue[10][5]=1.358;
+chartvalue[10][6]=1.937;
+ chartvalue[10][7]=1.798;
+ chartvalue[11][1]=0.347; chartvalue[11][2]=0.275;
+chartvalue[11][3]=0.679;
+ chartvalue[11][4]=0.395; chartvalue[11][5]=0.294;
+chartvalue[11][6]=0.579;
+ chartvalue[11][7]=0.213;
+ chartvalue[12][1]=0.262; chartvalue[12][2]=3.496;
+chartvalue[12][3]=3.108;
+ chartvalue[12][4]=0.998; chartvalue[12][5]=5.685;
+chartvalue[12][6]=2.494;
+ chartvalue[12][7]=3.048;
+ chartvalue[13][1]=0.03; chartvalue[13][2]=2.352;
+chartvalue[13][3]=2.268;
+ chartvalue[13][4]=0.237; chartvalue[13][5]=0.663;
+chartvalue[13][6]=1.62;
+ chartvalue[13][7]=1.448;
+ chartvalue[14][1]=0.179; chartvalue[14][2]=2.114;
+chartvalue[14][3]=1.778;
+ chartvalue[14][4]=0.631; chartvalue[14][5]=2.55;
+chartvalue[14][6]=1.578;
+ chartvalue[14][7]=2.526;
+ chartvalue[15][1]=0.835; chartvalue[15][2]=1.475;
+chartvalue[15][3]=1.534;
+ chartvalue[15][4]=0.039; chartvalue[15][5]=1.722;
+chartvalue[15][6]=2.456;
+ chartvalue[15][7]=2.280;
+ chartvalue[16][1]=0.382; chartvalue[16][2]=0.583;
+chartvalue[16][3]=1.052;
+ chartvalue[16][4]=0.419; chartvalue[16][5]=0.525;
+chartvalue[16][6]=0.916;
+ chartvalue[16][7]=0.628;
+ chartvalue[17][1]=0.169; chartvalue[17][2]=0.702;
+chartvalue[17][3]=0.955;
+ chartvalue[17][4]=0.654; chartvalue[17][5]=0.791;
+chartvalue[17][6]=0.843;
+ chartvalue[17][7]=0.647;
+ chartvalue[18][1]=0.824; chartvalue[18][2]=0.022;
+chartvalue[18][3]=0.308;
+ chartvalue[18][4]=0.152; chartvalue[18][5]=0.180;
+chartvalue[18][6]=0.156;
+ chartvalue[18][7]=0.044;
+ chartvalue[19][1]=0.24; chartvalue[19][2]=0;
+chartvalue[19][3]=0.00;
+ chartvalue[19][4]=0.456; chartvalue[19][5]=0.019;
+chartvalue[19][6]=0.00;
+ chartvalue[19][7]=0.00;
+ chartvalue[20][1]=0.00; chartvalue[20][2]=0.008; chartvalue[20][3]=0;
+ chartvalue[20][4]=0.013; chartvalue[20][5]=0.0; chartvalue[20][6]=0;
+ chartvalue[20][7]=0;
+}
+
+
+void Get_info_old()
+{
+ char junk;
+
+ printf("COILS version 2.1\n");
+ printf("ENTER INPUT FILE: ");
+ fflush(stdout);
+ scanf("%s",&inputfile[0]);
+ printf("-->%s\n",inputfile);
+ protein=fopen(inputfile,"r"); /* modify for reading */
+ if (protein==NULL) {
+ perror("problem");
+ printf("'%s' can not be found. Try again after verifying the file name.\n",inputfile);
+ exit(0);
+ } /* endif */
+
+ printf("ENTER OUTPUT FILE: ");
+ scanf("%s",&outputfile[0]);
+ pickup=fopen(outputfile,"w");
+ if (pickup==NULL) {
+ printf("The system could not open the output file\n");
+ printf(". Check the free space on disk\n");
+ exit(0);
+ } /* endif */
+ printf("Two scoring matrices are available:\n");
+ printf(" 1 - MTK\n");
+ printf(" 2 - MTIDK\n");
+ printf("Which matix? (1|2) <enter>? ");
+ chartchar='\0';
+ while ((chartchar!='1') && (chartchar!='2') && (chartchar!=(char) 13))
+ scanf("%c",&chartchar);
+ printf("Do you want a weight of 2.5 for positions a & d? (Y|N) <enter>? ");
+ weightchar='\0';
+ while ((weightchar!='n') && (weightchar!='y') && (weightchar!='Y') &&
+(weightchar!='N') && (weightchar!=(char) 13))
+ scanf("%c",&weightchar);
+
+ if(verbose) fprintf(pickup,"COILS version 2.1\n");
+ if (chartchar=='2')
+ {
+ if(verbose) fprintf(pickup,"using MTIDK matrix.\n");
+ Init_newchart();
+ } else {
+ if(verbose) fprintf(pickup,"using MTK matrix.\n");
+ Init_oldchart();
+ }
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ if(verbose) fprintf(pickup,"weights: a,d=2.5 and b,c,e,f,g=1.0\n");
+ hept_weight=10;
+ ad_weight=2.5;
+ } else {
+ if(verbose) fprintf(pickup,"no weights\n");
+ hept_weight=7;
+ ad_weight=1.0;
+ }
+ if(verbose) {
+ if(inputfile[0]!='-') fprintf(pickup,"Input file is %s\n",inputfile);
+ else fprintf(pickup,"Reading from standard input\n");
+ }
+}
+void Select_option()
+{
+ printf("OUTPUT OPTIONS:\n");
+ printf(" p - probabilities in columns\n");
+ printf(" a - probabilities in rows, abbreviated to the first digit\n");
+ printf(" b - scores (size of the scanning window defined by the user)\n");
+ printf(" c - only probabilities above a user-defined cutoff\n");
+ printf("ENTER OPTION: ");
+ while ((option!='a') && (option!='b') && (option!='c') &&
+(option!='A')
+&& (option!='B') && (option!='C') && (option!='p') && (option!='P')) {
+ scanf("%c",&option);
+ } /* endwhile */
+ printf("\n");
+ if ((option=='b') || (option=='B'))
+ {
+ printf("ENTER WINDOW SIZE: ");
+ scanf("%d",window);
+ } /* end */
+ else window=14;
+}
+
+
+int checkforend()
+{
+ int junklen;
+ char *last4;
+ int checkforend=0;
+
+ junkline[strlen(junkline)-1]='\0';
+ junklen=strlen(junkline);
+ if (junklen > 0)
+ {
+ if (junkline[0] == '>')
+ checkforend=1;
+ strcpy(nameprot,junkline);
+ }
+ if (junklen > 3)
+ {
+ last4=(char *) &(junkline[junklen-4]);
+ if (strcmp(last4," ..")==0)
+ {
+ checkforend=1;
+ strcpy(nameprot,junkline);
+ }
+ }
+ return checkforend;
+}
+void Skip_opening_info()
+ /* Skip opening lines of protein to get to amino acid sequence */
+{
+ int endingfound;
+
+ endingfound=checkforend();
+ while (!endingfound)
+ {
+ if (!feof(protein))
+ {
+ fgets(junkline,9999,protein);
+ endingfound=checkforend();
+ } else {
+ strcpy(nameprot," ");
+ fclose(protein);
+ protein=fopen(inputfile,"r");
+ endingfound=1;
+ } /* endif */
+ } /* endwhile */
+}
+void Store_residues()
+ /* stores only the amino acids of the sequence, skipping numbers,
+blanks,
+ and carriage returns and ending with '/' || '*' */
+{
+ char tempstore;
+ int validaa;
+
+ res_total=0;
+ tempstore=fgetc(protein);
+ while ((tempstore != '/') && (tempstore != '*') && (!feof(protein)))
+ {
+ switch (tempstore & 95) {
+ case 'L':
+ case 'I':
+ case 'V':
+ case 'M':
+ case 'F':
+ case 'Y':
+ case 'G':
+ case 'A':
+ case 'K':
+ case 'P':
+ case 'R':
+ case 'H':
+ case 'E':
+ case 'D':
+ case 'Q':
+ case 'N':
+ case 'S':
+ case 'T':
+ case 'C':
+ case 'W':
+ validaa=1;
+ break;
+ default:
+ validaa=0;
+ break;
+ } /* endswitch */
+
+ if (validaa)
+ {
+ res_total=res_total+1;
+ residue[res_total]=tempstore;
+ } /* endif */
+ tempstore=fgetc(protein);
+ } /* endwhile */
+ if (tempstore == '/') fgets(junkline,9999,protein); /* read other '/'
+at end */
+}
+
+
+int Getx(char res)
+{
+ int X;
+
+ switch (res & (char) 95) {
+ case 'L':
+ X=1; break;
+ case 'I':
+ X=2; break;
+ case 'V':
+ X=3; break;
+ case 'M':
+ X=4; break;
+ case 'F':
+ X=5; break;
+ case 'Y':
+ X=6; break;
+ case 'G':
+ X=7; break;
+ case 'A':
+ X=8; break;
+ case 'K':
+ X=9; break;
+ case 'R':
+ X=10; break;
+ case 'H':
+ X=11; break;
+ case 'E':
+ X=12; break;
+ case 'D':
+ X=13; break;
+ case 'Q':
+ X=14; break;
+ case 'N':
+ X=15; break;
+ case 'S':
+ X=16; break;
+ case 'T':
+ X=17; break;
+ case 'C':
+ X=18; break;
+ case 'W':
+ X=19; break;
+ case 'P':
+ X=20; break;
+ default:
+ X=-1; break;
+ } /* endswitch */
+ return X;
+}
+void Calculate(int startwindow,int endindex)
+ /* calculates best scores for each window frame */
+{
+ int X,x,y,extras,heptad_pos,window_pos,res_pos;
+ double root_inverse,scores,misc;
+ int hept,index,window;
+ double weight;
+ int root;
+
+ if(verbose) printf("Calculating...\n");
+ index=1;
+ window=startwindow;
+ while (index<=endindex)
+ {
+ if(verbose) printf(".");
+ fflush(stdout);
+ res_pos=0;
+ root=(window / 7) * hept_weight;
+ do {
+ res_pos=res_pos+1;
+ tempcalc[index][res_pos]=0;
+ for (heptad_pos=0; heptad_pos<=6; heptad_pos++) /* go through each residue in each heptad pos */
+ {
+ scores=1.0;
+ hept=1;
+ for (window_pos=0;window_pos<=window-1;window_pos++) /* get values at all 21 positions */
+ {
+ switch (residue[window_pos + res_pos] & 95){
+ case 'L':
+ x=1; break;
+ case 'I':
+ x=2; break;
+ case 'V':
+ x=3; break;
+ case 'M':
+ x=4; break;
+ case 'F':
+ x=5; break;
+ case 'Y':
+ x=6; break;
+ case 'G':
+ x=7; break;
+ case 'A':
+ x=8; break;
+ case 'K':
+ x=9; break;
+ case 'R':
+ x=10; break;
+ case 'H':
+ x=11; break;
+ case 'E':
+ x=12; break;
+ case 'D':
+ x=13; break;
+ case 'Q':
+ x=14; break;
+ case 'N':
+ x=15; break;
+ case 'S':
+ x=16; break;
+ case 'T':
+ x=17; break;
+ case 'C':
+ x=18; break;
+ case 'W':
+ x=19; break;
+ case 'P':
+ x=20; break;
+ default:
+ x=-1; break;
+ } /* endswitch */
+ y=(int) fmod((double)(window_pos + res_pos + heptad_pos),7.0);
+ if (y==0) y=7;
+ if (window_pos==0) hept=y;
+ if (y==0) y=7;
+ if ((y==1) || (y==4)) weight=ad_weight;
+ else weight=1.0;
+ root_inverse=1.0/(double) root;
+ misc=pow(chartvalue[x][y],weight);
+ scores=scores*(pow(misc,root_inverse));
+ } /* end of window_pos loop */
+ if (scores>tempcalc[index][res_pos])
+ {
+ tempcalc[index][res_pos]=scores;
+ temphept[index][res_pos]=(int) fmod((double) (hept-1),7.0);
+ }
+ } /* end of heptad_pos loop*/
+ } while (res_pos+window != res_total+1);
+ for (extras=1; extras<=window-1;extras++)
+ {
+ tempcalc[index][res_pos+extras]=tempcalc[index][res_pos];
+ temphept[index][res_pos+extras]=(int) fmod((double)
+(temphept[index][res_pos]+extras),7.0);
+ }
+ window=window+7;
+ index=index+1;
+ } /* for window sizes 14, 21 & 28 */
+ /*maximize loop*/
+ index=1;
+ window=startwindow;
+ while (index<=endindex)
+ {
+ res_pos=0;
+ do {
+ res_pos=res_pos+1;
+ calcnumb[index][res_pos]=tempcalc[index][res_pos];
+ heptnum[index][res_pos]=temphept[index][res_pos];
+ window_pos=0;
+ do {
+ window_pos=window_pos+1;
+ if (res_pos-window_pos<1) window_pos=window-1;
+ else
+ if (tempcalc[index][res_pos-window_pos]>calcnumb[index][res_pos])
+ {
+ calcnumb[index][res_pos]=tempcalc[index][res_pos-window_pos];
+ heptnum[index][res_pos]=(int)fmod((temphept[index][res_pos-window_pos]+window_pos),7.0);
+ } /* endif */
+ } while (window_pos!=window-1); /* enddo */
+ } while (res_pos!=res_total); /* enddo */
+ index=index+1;
+ window=window+7;
+ }
+ if(verbose) printf("\n");
+}
+double Calcprob(double x,double meancc,double stddevcc,double
+meangl,double stddevgl,double ratio_gl_cc)
+{
+ double prob1,prob2,prob3,prob4;
+
+ prob1=(0.5) * pow(((x-meancc) / stddevcc),2.0);
+ prob2=(0.5) * pow(((x-meangl) / stddevgl),2.0);
+ prob3=stddevgl * exp(-prob1);
+ prob4=ratio_gl_cc * stddevcc * exp(-prob2);
+ return (prob3) / (prob3+prob4);
+}
+char frame(int heptnum)
+{
+ switch (heptnum) {
+ case 0:
+ return 'a'; break;
+ case 1:
+ return 'b'; break;
+ case 2:
+ return 'c'; break;
+ case 3:
+ return 'd'; break;
+ case 4:
+ return 'e'; break;
+ case 5:
+ return 'f'; break;
+ case 6:
+ return 'g'; break;
+ } /* endswitch */
+ return 'x';
+}
+int i_trunc(double n)
+{
+ return (int) floor(n);
+}
+void Column_probs()
+{
+ int res_pos;
+ char fr14,fr21,fr28;
+ double prob14,prob21,prob28,old14,old21,old28,xx=1,lastxx=1;
+ int prevres,color=8,c14=9,c21=15,c28=12;
+ char t1[]="14",t2[]="21",t3[]="28";
+
+ fprintf(pickup,"%s\n",nameprot);
+ fprintf(pickup," Residue Window=14 Window=21 Window=28\n");
+ fprintf(pickup," Score Probability Score Probability Score Probability\n");
+ prevres=1;
+ for (res_pos=1;res_pos<=res_total;res_pos++)
+ {
+ fr14=frame(heptnum[1][res_pos]);
+ fr21=frame(heptnum[2][res_pos]);
+ fr28=frame(heptnum[3][res_pos]);
+ if (chartchar=='2')
+ {
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ prob28=Calcprob(calcnumb[3][res_pos],1.74,0.20,0.86,0.18,30);
+ prob21=Calcprob(calcnumb[2][res_pos],1.79,0.24,0.92,0.22,25);
+ prob14=Calcprob(calcnumb[1][res_pos],1.89,0.30,1.04,0.27,20);
+ } else {
+ prob28=Calcprob(calcnumb[3][res_pos],1.69,0.18,0.80,0.18,30);
+ prob21=Calcprob(calcnumb[2][res_pos],1.74,0.23,0.86,0.21,25);
+ prob14=Calcprob(calcnumb[1][res_pos],1.82,0.28,0.95,0.26,20);
+ } /* end */
+ } /* end */
+ else
+ {
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ if (calcnumb[3][res_pos] > 0.0)
+ prob28=Calcprob(calcnumb[3][res_pos],1.70,0.24,0.79,0.23,30);
+ else {
+ prob28=0.0;
+ fr28='x';
+ } /* endif */
+ if (calcnumb[2][res_pos] > 0.0)
+ prob21=Calcprob(calcnumb[2][res_pos],1.76,0.28,0.86,0.26,25);
+ else { prob21=0.0; fr21='x'; }
+ if (calcnumb[1][res_pos] > 0.0)
+ prob14=Calcprob(calcnumb[1][res_pos],1.88,0.34,1.00,0.33,20);
+ else { prob14=0.0; fr14='x'; }
+ } /* end */
+ else
+ {
+ if (calcnumb[3][res_pos] > 0.0)
+
+prob28=Calcprob(calcnumb[3][res_pos],1.628,0.243,0.770,0.202,30);
+ else { prob28=0.0; fr28='x'; }
+ if (calcnumb[2][res_pos] > 0.0)
+
+prob21=Calcprob(calcnumb[2][res_pos],1.683,0.285,0.828,0.236,25);
+ else { prob21=0.0; fr21='x'; }
+ if (calcnumb[1][res_pos] > 0.0)
+
+prob14=Calcprob(calcnumb[1][res_pos],1.782,0.328,0.936,0.289,20);
+ else { prob14=0.0; fr14='x'; }
+ } /* end */
+ }
+ if ((prob14>=peakmin) || (prob21>=peakmin) || (prob28>=peakmin))
+ {
+ if (res_pos-1!=prevres) fprintf(pickup,"...\n");
+ prevres=res_pos;
+ fprintf(pickup,"%5d %c",res_pos,residue[res_pos]);
+ fprintf(pickup," %c %5.3lf %5.3lf",fr14,calcnumb[1][res_pos],prob14);
+ fprintf(pickup," %c %5.3lf %5.3lf",fr21,calcnumb[2][res_pos],prob21);
+ fprintf(pickup," %c %5.3lf %5.3lf",fr28,calcnumb[3][res_pos],prob28);
+ fprintf(pickup,"\n");
+ }
+ }
+}
+void Row_probs()
+{
+ int res_pos,maxline,startpos;
+
+ fprintf(pickup,"%s\n",nameprot);
+ for (res_pos=1;res_pos<=res_total;res_pos++)
+ {
+ fr14[res_pos]=frame(heptnum[1][res_pos]);
+ fr21[res_pos]=frame(heptnum[2][res_pos]);
+ fr28[res_pos]=frame(heptnum[3][res_pos]);
+ if (chartchar=='2')
+ {
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ pr14[res_pos]=i_trunc(Calcprob(calcnumb[1][res_pos],1.89,0.30,1.04,0.27,20)*10.0);
+ pr21[res_pos]=i_trunc(Calcprob(calcnumb[2][res_pos],1.79,0.24,0.92,0.22,25)*10.0);
+ pr28[res_pos]=i_trunc(Calcprob(calcnumb[3][res_pos],1.74,0.20,0.86,0.18,30)*10.0);
+ } /* end */
+ else
+ {
+ pr14[res_pos]=i_trunc(Calcprob(calcnumb[1][res_pos],1.82,0.28,0.95,0.26,20)*10.0);
+ pr21[res_pos]=i_trunc(Calcprob(calcnumb[2][res_pos],1.74,0.23,0.86,0.21,25)*10.0);
+ pr28[res_pos]=i_trunc(Calcprob(calcnumb[3][res_pos],1.69,0.18,0.80,0.18,30)*10.0);
+ } /* end */
+ } /* end */
+ else
+ {
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ if (calcnumb[3][res_pos] > 0)
+ pr28[res_pos]=i_trunc(Calcprob(calcnumb[3][res_pos],1.70,0.24,0.79,0.23,30)*10.0);
+ else { pr28[res_pos]=0; fr28[res_pos]='x'; }
+ if (calcnumb[2][res_pos] > 0)
+ pr21[res_pos]=i_trunc(Calcprob(calcnumb[2][res_pos],1.76,0.28,0.86,0.26,25)*10.0);
+ else { pr21[res_pos]=0; fr21[res_pos]='x'; }
+ if (calcnumb[1][res_pos] > 0)
+ pr14[res_pos]=i_trunc(Calcprob(calcnumb[1][res_pos],1.88,0.34,1.00,0.33,20)*10.0);
+ else { pr14[res_pos]=0; fr14[res_pos]='x'; }
+ } /* end */
+ else
+ {
+ if (calcnumb[3][res_pos] > 0)
+ pr28[res_pos]=i_trunc(Calcprob(calcnumb[3][res_pos],1.628,0.243,0.770,0.202,30)*10.0);
+ else { pr28[res_pos]=0; fr28[res_pos]='x'; }
+ if (calcnumb[2][res_pos] > 0)
+ pr21[res_pos]=i_trunc(Calcprob(calcnumb[2][res_pos],1.683,0.285,0.828,0.236,25)*10.0);
+ else { pr21[res_pos]=0; fr21[res_pos]='x'; }
+ if (calcnumb[1][res_pos] > 0)
+ pr14[res_pos]=i_trunc(Calcprob(calcnumb[1][res_pos],1.782,0.328,0.936,0.289,20)*10.0);
+ else { pr14[res_pos]=0; fr14[res_pos]='x'; }
+ } /* end */
+ }
+ if (pr14[res_pos]>9) pr14[res_pos]=9;
+ if (pr21[res_pos]>9) pr21[res_pos]=9;
+ if (pr28[res_pos]>9) pr28[res_pos]=9;
+ }
+ maxline=60;
+ startpos=1;
+ while (startpos+maxline-1 <= res_total)
+ {
+ fprintf(pickup,"%d\n",startpos);
+ fprintf(pickup," . | . | . | . | . | . |\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%c",residue[res_pos]);
+ fprintf(pickup,"\n\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%d",pr14[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%d",pr21[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%d",pr28[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%c",fr14[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%c",fr21[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=startpos+maxline-1;res_pos++)
+ fprintf(pickup,"%c",fr28[res_pos]);
+ fprintf(pickup,"\n");
+ fprintf(pickup,"\n");
+ startpos=startpos+60;
+ }
+ /* finish up! */
+ if (startpos<=res_total)
+ {
+ fprintf(pickup,"%d\n",startpos);
+ fprintf(pickup," . | . | . | . | . | . |\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%c",residue[res_pos]);
+ fprintf(pickup,"\n");
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%d",pr14[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%d",pr21[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%d",pr28[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%c",fr14[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%c",fr21[res_pos]);
+ fprintf(pickup,"\n");
+ for (res_pos=startpos;res_pos<=res_total;res_pos++)
+ fprintf(pickup,"%c",fr28[res_pos]);
+ fprintf(pickup,"\n");
+ fprintf(pickup,"\n");
+ }
+}
+
+
+void Scores_only()
+{
+ int res_pos;
+ char fr;
+
+ fprintf(pickup,"%s\n",nameprot);
+ fprintf(pickup," Residue Frame Score\n");
+ for (res_pos=1;res_pos<=res_total;res_pos++)
+ {
+ if (calcnumb[1][res_pos] > 0.0)
+ {
+ fr=frame(heptnum[1][res_pos]);
+ fprintf(pickup," %6d %c %c %2.6lf\n",res_pos,residue[res_pos],fr, calcnumb[1][res_pos]);
+ } else
+ fprintf(pickup," %6d %c x %2.6lf\n",res_pos,residue[res_pos], calcnumb[1][res_pos]);
+ }
+}
+
+void Calc_print()
+{
+ int num_wind;
+
+ if (res_total >= window)
+ {
+ num_wind=3;
+ if ((option=='a') || (option=='A'))
+ {
+ Calculate(window,num_wind);
+ Row_probs();
+ } else
+ if ((option=='b') || (option=='B'))
+ {
+ num_wind=1;
+ Calculate(window,num_wind);
+ Scores_only();
+ } /* end */
+ else
+ {
+ Calculate(window,num_wind);
+ Column_probs();
+ }
+ } else
+ fprintf(pickup,"Protein too short. Length=%d\n",res_total);
+}
+
+void exit_error() {
+ fprintf(stderr,"format: coils [options] < [fasta file] > [output file]\n");
+ fprintf(stderr," -mtk (use MTK matrix) \n");
+ fprintf(stderr," -mtidk (use MTIDK matrix) \n");
+ fprintf(stderr," -ad2.5 (a weight of 2.5 for positions a & d)\n");
+ fprintf(stderr," -w <window size> [default is 21]\n");
+ fprintf(stderr," -pm <minimum probabilty> [default is 0.0]\n");
+ fprintf(stderr," -v (verbose mode)\n");
+ fprintf(stderr," -out <char>\n");
+ fprintf(stderr," p = propabilities in columns \n");
+ fprintf(stderr," a = probabilities in rows\n");
+ fprintf(stderr," b = scores (size of the scanning window defined by the user)\n");
+ exit(-1);
+}
+
+void Get_info(int argc, char *argv[])
+{
+
+ /* Modified version of Get_info()
+ * Rob Russell, March 1999 makes commands command line driven */
+
+ /* Defaults */
+ peakmin = 0.0;
+ window = 14;
+ option = 'p';
+ verbose = 0;
+ hept_weight =7;
+ ad_weight = 1.0;
+
+ chartchar = '1';
+ Init_oldchart();
+
+ strcpy(&inputfile[0],"-");
+ for(i=1; i<argc; ++i) {
+ if(argv[i][0]!='-') {
+ strcpy(&inputfile[0],argv[i]);
+ } else if(strcmp(&argv[i][1],"v")==0) {
+ verbose = 1;
+ } else if(strcmp(&argv[i][1],"mtk")==0) {
+ chartchar = '1';
+ Init_oldchart();
+ } else if(strcmp(&argv[i][1],"mtidk")==0) {
+ chartchar = '2';
+ Init_newchart();
+ } else if(strcmp(&argv[i][1],"ad2.5")==0) {
+ weightchar = 'y';
+ hept_weight=10;
+ ad_weight=2.5;
+ } else if(strcmp(&argv[i][1],"pm")==0) {
+ if((i+1)>=argc) { exit_error(); }
+ sscanf(argv[i+1],"%f",&peakmin);
+ i++;
+ } else if(strcmp(&argv[i][1],"w")==0) {
+ if((i+1)>=argc) { exit_error(); }
+ sscanf(argv[i+1],"%d",&window);
+ i++;
+ } else if(strcmp(&argv[i][1],"out")==0) {
+ if((i+1)>=argc) { exit_error(); }
+ option = argv[i+1][0];
+ if((option!='p') && (option!='a') && (option!='b') && (option!='c')) {
+ exit_error();
+ }
+ i++;
+ } else {
+ exit_error();
+ }
+ }
+
+ if(inputfile[0]!='-') {
+ protein=fopen(inputfile,"r"); /* modify for reading */
+ if (protein==NULL) {
+ perror("problem");
+ fprintf(stderr,"Error: file %s does not exist\n",inputfile);
+ exit(0);
+ }
+ if(verbose) printf("Input file is: %s\n",inputfile);
+ } else {
+ protein = stdin;
+ if(verbose) printf("Reading from standard input\n");
+ }
+ pickup = stdout;
+ if(verbose) fprintf(pickup,"COILS version 2.1\n");
+ if (chartchar=='2')
+ {
+ if(verbose) fprintf(pickup,"using MTIDK matrix.\n");
+ Init_newchart();
+ } else {
+ if(verbose) fprintf(pickup,"using MTK matrix.\n");
+ Init_oldchart();
+ }
+ if ((weightchar=='y') || (weightchar=='Y'))
+ {
+ if(verbose) fprintf(pickup,"weights: a,d=2.5 and b,c,e,f,g=1.0\n");
+ hept_weight=10;
+ ad_weight=2.5;
+ } else {
+ if(verbose) fprintf(pickup,"no weights\n");
+ hept_weight=7;
+ ad_weight=1.0;
+ }
+}
+
+
+main(int argc, char *argv[], char *envp[])
+{ /* Main Program */
+ Get_info(argc, argv);
+
+/* Get_info_old();
+ Select_option();
+*/
+ fgets(junkline,9999,protein);
+
+ do {
+ Skip_opening_info();
+ if (!nomore)
+ {
+ Store_residues();
+ Calc_print();
+ if (!feof(protein)) fgets(junkline,9999,protein);
+ }
+ } while (!feof(protein)); /* enddo */
+ fclose(protein);
+ fclose(pickup);
+}
+