1 /* $Id: res_stats.c 1227 2013-09-26 19:19:28Z wrp $ */
2
3 /* copyright (C) 2005, 2014 by William R. Pearson and The Rector &
4 Vistors of the University of Virginia */
5
6 /* Licensed under the Apache License, Version 2.0 (the "License");
7 you may not use this file except in compliance with the License.
8 You may obtain a copy of the License at
9
10 http://www.apache.org/licenses/LICENSE-2.0
11
12 Unless required by applicable law or agreed to in writing,
13 software distributed under this License is distributed on an "AS
14 IS" BASIS, WITHOUT WRRANTIES OR CONDITIONS OF ANY KIND, either
15 express or implied. See the License for the specific language
16 governing permissions and limitations under the License.
17 */
18
19 /* calculate stats from results file using scalesws.c */
20
21 #include <stdio.h>
22 #include <stdlib.h>
23 #include <string.h>
24
25 #include <limits.h>
26 #include <math.h>
27
28 #define MAX_LLEN 200
29
30 #define LN_FACT 10.0
31
32 #include "defs.h"
33 #include "structs.h"
34 #include "param.h"
35
36 struct beststr {
37 int score; /* smith-waterman score */
38 int sscore; /* duplicate for compatibility with fasta */
39 double comp;
40 double H;
41 double zscore;
42 double escore;
43 int n1;
44 #ifndef USE_FTELLO
45 long lseek; /* position in library file */
46 #else
47 off_t lseek;
48 #endif
49 int cont; /* offset into sequence */
50 int frame;
51 int lib;
52 char libstr[13];
53 } *bbp, *bestptr, **bptr, *best;
54
55 struct stat_str {
56 int score;
57 int n1;
58 double comp;
59 double H;
60 };
61
62 static struct db_str qtt = {0l, 0l, 0};
63
64 char info_gstring2[MAX_STR]; /* string for label */
65 char info_gstring3[MAX_STR];
66 char info_hstring1[MAX_STR];
67
68 FILE *outfd;
69
70 int nbest; /* number of sequences better than bestcut in best */
71 int bestcut=1; /* cut off for getting into MAX_BEST */
72 int bestfull;
73
74 int dohist = 0;
75 int zsflag = 1;
76 int outtty=1;
77 int llen=40;
78
79 /* statistics functions */
80 extern void
81 process_hist(struct stat_str *sptr, int nstat, struct pstruct pst,
82 struct hist_str *hist, void **);
83 extern void addhistz(double, struct hist_str *); /* scaleswn.c */
84 void selectbestz(struct beststr **, int, int );
85
86 extern double zs_to_E(double, int, int, long, struct db_str);
87 extern double zs_to_Ec(double zs, long entries);
88
89 extern double find_z(int score, int length, double comp, void *);
90
91 void prhist(FILE *, struct mngmsg, struct pstruct, struct hist_str,
92 int, struct db_str, char *);
93
94 int nshow=20, mshow=50, ashow= -1;
95 double e_cut=10.0;
96
main(argc,argv)97 main(argc, argv)
98 int argc; char **argv;
99 {
100 FILE *fin;
101 char line[512];
102 int max, icol, iarg, i, qsfnum, lsfnum, n0, n1, s[3], frame;
103 double comp, H;
104 int idup, ndup, max_s;
105 char libstr[MAX_UID], *bp;
106 char bin_file[80];
107 FILE *bout=NULL;
108 struct mngmsg m_msg; /* Message from host to manager */
109 struct pstruct pst;
110 struct stat_str *stats;
111 int nstats;
112 double zscor, mu, var;
113
114 #if defined(UNIX)
115 outtty = isatty(1);
116 #else
117 outtty = 1;
118 #endif
119
120 if (argc < 2 ) {
121 fprintf(stderr," useage - res_stats -c col -r bin_file file\n");
122 exit(1);
123 }
124
125 m_msg.db.length = qtt.length = 0l;
126 m_msg.db.entries = m_msg.db.carry = qtt.entries = qtt.carry = 0;
127 m_msg.pstat_void = NULL;
128 m_msg.hist.hist_a = NULL;
129 m_msg.nohist = 0;
130 m_msg.markx = 0;
131
132 pst.n0 = 200; /* sensible dummy value */
133 pst.zsflag = 1;
134 pst.dnaseq = 0;
135 pst.histint = 2;
136
137 bin_file[0]='\0';
138 icol = 1;
139 iarg = 1;
140 ndup = 1;
141 while (1) {
142 if (argv[iarg][0]=='-' && argv[iarg][1]=='c') {
143 sscanf(argv[iarg+1],"%d",&icol);
144 iarg += 2;
145 }
146 else if (argv[iarg][0]=='-' && argv[iarg][1]=='r') {
147 strncpy(bin_file,argv[iarg+1],sizeof(bin_file));
148 iarg += 2;
149 }
150 else if (argv[iarg][0]=='-' && argv[iarg][1]=='z') {
151 sscanf(argv[iarg+1],"%d",&pst.zsflag);
152 iarg += 2;
153 }
154 else if (argv[iarg][0]=='-' && argv[iarg][1]=='n') {
155 pst.dnaseq = 1;
156 iarg += 1;
157 }
158 else if (argv[iarg][0]=='-' && argv[iarg][1]=='s') {
159 sscanf(argv[iarg+1],"%d",&ndup);
160 iarg += 2;
161 }
162 else if (argv[iarg][0]=='-' && argv[iarg][1]=='q') {
163 outtty = 0;
164 iarg += 1;
165 }
166 else break;
167 }
168
169 icol--;
170
171 if ((fin=fopen(argv[iarg],"r"))==NULL) {
172 fprintf(stderr," cannot open %s\n",argv[1]);
173 exit(1);
174 }
175
176 if (bin_file[0]!='\0' && ((bout=fopen(bin_file,"w"))==NULL)) {
177 fprintf(stderr,"cannot open %s for output\n",bin_file);
178 }
179
180 if ((stats =
181 (struct stat_str *)malloc((MAX_STATS)*sizeof(struct stat_str)))==NULL)
182 s_abort ("Cannot allocate stats struct","");
183 nstats = 0;
184
185 initbest(MAX_BEST+1); /* +1 required for select() */
186
187 for (nbest=0; nbest<MAX_BEST+1; nbest++)
188 bptr[nbest] = &best[nbest];
189 bptr++; best++;
190 best[-1].score= BIGNUM;
191
192 nbest = 0;
193
194 pst.Lambda=0.232;
195 pst.K = 0.11;
196 pst.H = 0.34;
197
198 /* read the best scores from the results file */
199
200 max_s = -1;
201 idup = 0;
202
203 /* get first line with sequence length */
204 fgets(line,sizeof(line),fin);
205 sscanf(line,"%d",&n0);
206 if (n0 > 0) pst.n0 = n0;
207
208 while (fgets(line,sizeof(line),fin)!=NULL) {
209 if (line[0]=='/' && line[1]=='*') {
210 fputs(line,stdout);
211 strncpy(info_gstring2,line,sizeof(info_gstring2));
212 if ((bp=strchr(info_gstring2,'\n'))!=NULL) *bp = '\0';
213 break;
214 }
215 if (line[0]==';') {
216 if ((bp=strchr(line,'|'))!=NULL) qsfnum = atoi(bp+1);
217 else continue;
218 if ((bp=strchr(line,'('))!=NULL) {
219 n0 = atoi(bp+1);
220 pst.n0 = n0;
221 }
222 else {
223 fprintf(stderr, "cannot find n0:\n %s\n",line);
224 continue;
225 }
226 }
227 else {
228 sscanf(line,"%s %d %d %d %lf %lf %d %d %d",
229 libstr,&lsfnum,&n1,&frame,&comp, &H, &s[0],&s[1],&s[2]);
230 if (lsfnum==0 && n1==0) {
231 fputs(line,stderr);
232 continue;
233 }
234 if (n1 < 10 || s[icol]<=0) fputs(line,stderr);
235 idup++;
236
237 if (s[icol] > max_s) max_s = s[icol];
238 if (idup < ndup) continue;
239
240 m_msg.db.entries++;
241 m_msg.db.length += n1;
242
243 if (dohist) addhistz(zscor=find_z(max_s,n1,comp,m_msg.pstat_void),
244 &m_msg.hist);
245 else zscor = (double)max_s;
246
247 if (nstats < MAX_STATS) {
248 stats[nstats].n1 = n1;
249 stats[nstats].comp = comp;
250 stats[nstats].H = H;
251 stats[nstats++].score = max_s;
252 }
253
254 else if (!dohist) {
255 /* do_bout(bout,stats,nstats); */
256 process_hist(stats,nstats,pst,&m_msg.hist, &m_msg.pstat_void);
257 for (i=0; i<nbest; i++)
258 bptr[i]->zscore =
259 find_z(bptr[i]->score,bptr[i]->n1,bptr[i]->comp,
260 m_msg.pstat_void);
261 dohist = 1;
262 }
263
264 if (dohist) {
265 zscor =find_z(max_s,n1,comp,m_msg.pstat_void);
266 addhistz(zscor,&m_msg.hist);
267 }
268 else zscor = (double)max_s;
269
270 if (nbest >= MAX_BEST) {
271 bestfull = nbest-MAX_BEST/4;
272 selectz(bestfull-1,nbest);
273 bestcut = (int)(bptr[bestfull-1]->zscore+0.5);
274 nbest = bestfull;
275 }
276 bestptr = bptr[nbest];
277 bestptr->score = max_s;
278 bestptr->sscore = max_s;
279 bestptr->n1 = n1;
280 bestptr->comp = comp;
281 bestptr->H = H;
282 bestptr->lib = lsfnum;
283 bestptr->zscore = zscor;
284 strncpy(bestptr->libstr,libstr,12);
285 bestptr->libstr[12]='\0';
286 nbest++;
287
288 max_s = -1;
289 idup = 0;
290 }
291 } /* done with reading results */
292
293 if (!dohist) {
294 if (nbest < 20) {
295 zsflag = 0;
296 }
297 else {
298 /* do_bout(bout,stats,nstats); */
299 process_hist(stats,nstats,pst,&m_msg.hist,&m_msg.pstat_void);
300 for (i=0; i<nbest; i++)
301 bptr[i]->zscore =
302 find_z(bptr[i]->score,bptr[i]->n1,bptr[i]->comp,m_msg.pstat_void);
303 dohist = 1;
304 }
305 }
306
307 printf(" using n0: %d\n",pst.n0);
308
309 /* print histogram, statistics */
310
311 m_msg.nbr_seq = m_msg.db.entries;
312 pst.zdb_size = m_msg.db.entries;
313 /* get_param(&pst, info_gstring2,info_gstring3); */
314
315 prhist(stdout,m_msg,pst,m_msg.hist,nstats,m_msg.db,info_gstring2);
316
317 if (!zsflag) sortbest();
318 else {
319 sortbestz(bptr,nbest);
320 for (i=0; i<nbest; i++)
321 bptr[i]->escore = zs_to_E(bptr[i]->zscore,bptr[i]->n1,pst.dnaseq,
322 pst.zdb_size, m_msg.db);
323 }
324
325 outfd = stdout;
326 showbest(m_msg.db); /* display best matches */
327 }
328
initbest(nbest)329 initbest(nbest) /* allocate arrays for best sort */
330 int nbest;
331 {
332
333 if ((best=(struct beststr *)calloc((size_t)nbest,sizeof(struct beststr)))
334 == NULL) {fprintf(stderr,"cannot allocate best struct\n"); exit(1);}
335 if ((bptr=(struct beststr **)calloc((size_t)nbest,sizeof(struct beststr *)))
336 == NULL) {fprintf(stderr,"cannot allocate bptr\n"); exit(1);}
337 }
338
339 void
prhist(FILE * fd,struct mngmsg m_msg,struct pstruct pst,struct hist_str hist,int nstats,struct db_str ntt,char * info_gstring2)340 prhist(FILE *fd, struct mngmsg m_msg,
341 struct pstruct pst,
342 struct hist_str hist,
343 int nstats,
344 struct db_str ntt,
345 char *info_gstring2)
346 {
347 int i,j,hl,hll, el, ell, ev;
348 char hline[80], pch, *bp;
349 int mh1, mht;
350 int maxval, maxvalt, dotsiz, ddotsiz,doinset;
351 double cur_e, prev_e, f_int;
352 double max_dev, x_tmp;
353 double db_tt;
354 int n_chi_sq, cum_hl, max_i;
355
356
357 fprintf(fd,"\n");
358
359 if (pst.zsflag < 0 || nstats <= 10) {
360 fprintf(fd, "%7ld residues in %5ld sequences\n", ntt.length,ntt.entries);
361 fprintf(fd,"\n%s\n",info_gstring2);
362 return;
363 }
364
365 max_dev = 0.0;
366 mh1 = hist.maxh-1;
367 mht = (3*hist.maxh-3)/4 - 1;
368
369 if (!m_msg.nohist && mh1 > 0) {
370 for (i=0,maxval=0,maxvalt=0; i<hist.maxh; i++) {
371 if (hist.hist_a[i] > maxval) maxval = hist.hist_a[i];
372 if (i >= mht && hist.hist_a[i]>maxvalt) maxvalt = hist.hist_a[i];
373 }
374 n_chi_sq = 0;
375 cum_hl = -hist.hist_a[0];
376 dotsiz = (maxval-1)/60+1;
377 ddotsiz = (maxvalt-1)/50+1;
378 doinset = (ddotsiz < dotsiz && dotsiz > 2);
379
380 if (pst.zsflag>=0)
381 fprintf(fd," opt E()\n");
382 else
383 fprintf(fd," opt\n");
384
385 prev_e = zs_to_Ec((double)(hist.min_hist-hist.histint/2),hist.entries);
386 for (i=0; i<=mh1; i++) {
387 pch = (i==mh1) ? '>' : ' ';
388 pch = (i==0) ? '<' : pch;
389 hll = hl = hist.hist_a[i];
390 if (pst.zsflag>=0) {
391 cum_hl += hl;
392 f_int = (double)(i*hist.histint+hist.min_hist)+(double)hist.histint/2.0;
393 cur_e = (double)zs_to_Ec(f_int,hist.entries);
394 ev = el = ell = (int)(cur_e - prev_e + 0.5);
395 if (hl > 0 && i > 5 && i < (90-hist.min_hist)/hist.histint) {
396 x_tmp = fabs(cum_hl - cur_e);
397 if ( x_tmp > max_dev) {
398 max_dev = x_tmp;
399 max_i = i;
400 }
401 n_chi_sq++;
402 }
403 if ((el=(el+dotsiz-1)/dotsiz) > 60) el = 60;
404 if ((ell=(ell+ddotsiz-1)/ddotsiz) > 40) ell = 40;
405 fprintf(fd,"%c%3d %5d %5d:",
406 pch,(i<mh1)?(i)*hist.histint+hist.min_hist :
407 mh1*hist.histint+hist.min_hist,hl,ev);
408 }
409 else fprintf(fd,"%c%3d %5d :",
410 pch,(i<mh1)?(i)*hist.histint+hist.min_hist :
411 mh1*hist.histint+hist.min_hist,hl);
412
413 if ((hl=(hl+dotsiz-1)/dotsiz) > 60) hl = 60;
414 if ((hll=(hll+ddotsiz-1)/ddotsiz) > 40) hll = 40;
415 for (j=0; j<hl; j++) hline[j]='=';
416 if (pst.zsflag>=0) {
417 if (el <= hl ) {
418 if (el > 0) hline[el-1]='*';
419 hline[hl]='\0';
420 }
421 else {
422 for (j = hl; j < el; j++) hline[j]=' ';
423 hline[el-1]='*';
424 hline[hl=el]='\0';
425 }
426 }
427 else hline[hl] = 0;
428 if (i==1) {
429 for (j=hl; j<10; j++) hline[j]=' ';
430 sprintf(&hline[10]," one = represents %d library sequences",dotsiz);
431 }
432 if (doinset && i == mht-2) {
433 for (j = hl; j < 10; j++) hline[j]=' ';
434 sprintf(&hline[10]," inset = represents %d library sequences",ddotsiz);
435 }
436 if (i >= mht&& doinset ) {
437 for (j = hl; j < 10; j++) hline[j]=' ';
438 hline[10]=':';
439 for (j = 11; j<11+hll; j++) hline[j]='=';
440 hline[11+hll]='\0';
441 if (pst.zsflag>=0) {
442 if (ell <= hll) hline[10+ell]='*';
443 else {
444 for (j = 11+hll; j < 10+ell; j++) hline[j]=' ';
445 hline[10+ell] = '*';
446 hline[11+ell] = '\0';
447 }
448 }
449 }
450
451 fprintf(fd,"%s\n",hline);
452 prev_e = cur_e;
453 }
454 }
455
456 if (ntt.carry==0) {
457 fprintf(fd, "%7ld residues in %5ld sequences\n", ntt.length, ntt.entries);
458 }
459 else {
460 db_tt = (double)ntt.carry*(double)LONG_MAX + (double)ntt.length;
461 fprintf(fd, "%.0f residues in %5ld library sequences\n", db_tt, ntt.entries);
462 }
463
464 if (pst.zsflag>=0) {
465 if (MAX_STATS < hist.entries)
466 fprintf(fd," statistics extrapolated from %d to %ld sequences\n",
467 MAX_STATS,hist.entries);
468 /* summ_stats(stat_info); */
469 fprintf(fd," %s\n",hist.stat_info);
470 if (!m_msg.nohist && cum_hl > 0)
471 fprintf(fd," Kolmogorov-Smirnov statistic: %6.4f (N=%d) at %3d\n",
472 max_dev/(double)cum_hl, n_chi_sq,max_i*hist.histint+hist.min_hist);
473 if (m_msg.markx & MX_M10FORM) {
474 while ((bp=strchr(hist.stat_info,'\n'))!=NULL) *bp=' ';
475 if (cum_hl <= 0) cum_hl = -1;
476 sprintf(info_hstring1,"; mp_extrap: %d %ld\n; mp_stats: %s\n; mp_KS: %6.4f (N=%d) at %3d\n",
477 MAX_STATS,hist.entries,hist.stat_info,max_dev/(double)cum_hl, n_chi_sq,max_i*hist.histint+hist.min_hist);
478 }
479 }
480 fprintf(fd,"\n%s\n",info_gstring2);
481 fflush(fd);
482 }
483
showbest(struct db_str ntt)484 showbest(struct db_str ntt)
485 {
486 int ib, istart, istop;
487 char bline[200], fmt[40], pad[200];
488 char rline[20];
489 int ntmp;
490 int lcont, ccont, loff;
491 int hcutoff;
492
493 sprintf(fmt,"%%-%ds (%%3d)",llen-10);
494
495 nshow = min(20,nbest);
496 mshow = min(20,nbest);
497
498 if (outtty) {
499 printf(" How many scores would you like to see? [%d] ",nshow);
500 fflush(stdout);
501 if (fgets(rline,sizeof(rline),stdin)==NULL) exit(0);
502 if (rline[0]!='\n' && rline[0]!=0) sscanf(rline,"%d",&nshow);
503 if (nshow<=0) nshow = min(20,nbest);
504 }
505 else nshow=mshow;
506
507 memset(pad,' ',llen-10);
508 pad[llen-31]='\0';
509 if (zsflag)
510 fprintf(outfd,"The best scores are:%s s-w Z-score E(%ld)\n",pad,ntt.entries);
511 else
512 fprintf(outfd,"The best scores are:%s s-w\n",pad);
513
514 if (outfd != stdout)
515 if (zsflag)
516 fprintf(stdout,"The best scores are:%s s-w Z-score E(%ld)\n",pad,ntt.entries);
517 else
518 fprintf(stdout,"The best scores are:%s s-w\n",pad);
519
520 istart = 0;
521 l1: istop = min(nbest,nshow);
522 for (ib=istart; ib<istop; ib++) {
523 bbp = bptr[ib];
524
525 if (!outtty && zsflag && bbp->escore > e_cut) {
526 nshow = ib;
527 goto done;
528 }
529
530 sprintf(bline,"%-12s %d",bbp->libstr,bbp->lib);
531 bline[13]='\0';
532
533 fprintf(outfd,fmt,bline,bbp->n1);
534
535 if (zsflag)
536 fprintf(outfd,"%4d %4.1f %6.2g\n",
537 bbp->score,bbp->zscore,
538 bbp->escore);
539 else
540 fprintf(outfd,"%4d\n",bbp->score);
541
542 if (outfd!=stdout) {
543 fprintf(stdout,fmt,bline,bbp->n1);
544 if (zsflag)
545 printf("%4d %4.1f %6.2g\n",
546 bbp->score,bbp->zscore,
547 bbp->escore);
548 else
549 printf("%4d\n",bbp->score);
550 }
551 }
552
553 fflush(outfd); if (outfd!=stdout) fflush(stdout);
554
555 if (outtty) {
556 printf(" More scores? [0] ");
557 fflush(stdout);
558 if (fgets(rline,sizeof(rline),stdin)==NULL) exit(0);
559 ntmp = 0;
560 if (rline[0]!='\n' && rline[0]!=0) sscanf(rline,"%d",&ntmp);
561 if (ntmp<=0) ntmp = 0;
562 if (ntmp>0) {
563 istart = istop;
564 nshow += ntmp;
565 mshow += ntmp;
566 goto l1;
567 }
568 }
569 else if (zsflag && bbp->escore < e_cut) {
570 istart=istop;
571 nshow += 10;
572 goto l1;
573 }
574
575 done:
576 if (outfd!=stdout) fprintf(outfd,"\n");
577 }
578
selectz(k,n)579 selectz(k,n) /* k is rank in array */
580 int k,n;
581 {
582 int t, i, j, l, r;
583 double v;
584 struct beststr *tmptr;
585
586 l=0; r=n-1;
587
588 while ( r > l ) {
589 i = l-1;
590 j = r;
591 v = bptr[r]->zscore;
592 do {
593 while (bptr[++i]->zscore > v ) ;
594 while (bptr[--j]->zscore < v ) ;
595 tmptr = bptr[i]; bptr[i]=bptr[j]; bptr[j]=tmptr;
596 } while (j > i);
597 bptr[j]=bptr[i]; bptr[i]=bptr[r]; bptr[r]=tmptr;
598 if (i>=k) r = i-1;
599 if (i<=k) l = i+1;
600 }
601 }
602
sortbest()603 sortbest()
604 {
605 int cmps(), cmp1(), cmpa(), cmpz();
606 ksort(bptr,nbest,cmps);
607 }
608
sortbeste()609 sortbeste()
610 {
611 int cmpe();
612 ksort(bptr,nbest,cmpe);
613 }
614
sortbestz()615 sortbestz()
616 {
617 int cmpz();
618 ksort(bptr,nbest,cmpz);
619 }
620
621 cmps(ptr1,ptr2)
622 struct beststr *ptr1, *ptr2;
623 {
624 if (ptr1->score < ptr2->score) return (1);
625 else if (ptr1->score > ptr2->score) return (-1);
626 else return (0);
627 }
628
629 cmpe(ptr1,ptr2)
630 struct beststr *ptr1, *ptr2;
631 {
632 if (ptr1->escore < ptr2->escore) return (-1);
633 else if (ptr1->escore > ptr2->escore) return (1);
634 else return (0);
635 }
636
637 cmpz(ptr1,ptr2)
638 struct beststr *ptr1, *ptr2;
639 {
640 if (ptr1->zscore < ptr2->zscore) return (1);
641 else if (ptr1->zscore > ptr2->zscore) return (-1);
642 else return (0);
643 }
644
ksort(v,n,comp)645 ksort(v,n,comp)
646 char *v[]; int n, (*comp)();
647 {
648 int gap, i, j;
649 char *tmp;
650
651 for (gap=n/2; gap>0; gap/=2)
652 for (i=gap; i<n; i++)
653 for (j=i-gap; j>=0; j -= gap) {
654 if ((*comp)(v[j],v[j+gap]) <=0)
655 break;
656 tmp = v[j]; v[j]=v[j+gap]; v[j+gap]=tmp;
657 }
658 }
659
660 /*
661 do_bout(FILE *bout,struct stat_str **bptr, int nbest)
662 {
663 int i, min_hist, max_hist;
664 double mu, var;
665
666 if (bout==NULL) return;
667
668 inithist();
669 for (i = 0; i<nbest; i++)
670 addhist(bptr[i]->score,bptr[i]->n1);
671
672 for (i=0; i<MAX_LLEN; i++)
673 if (llen_hist[i]>0) {
674 min_hist=i;
675 break;
676 }
677
678 for (i=MAX_LLEN-1; i>=0; i--)
679 if (llen_hist[i]>0) {
680 max_hist=i;
681 break;
682 }
683
684 for (i=min_hist; i<=max_hist; i++) {
685 mu=(double)score_sums[i]/(double)llen_hist[i];
686 if (llen_hist[i]>1) {
687 var = ((double)score2_sums[i]-(double)llen_hist[i]*mu*mu)/
688 (double)(llen_hist[i]-1);
689
690 fprintf(bout,"%d\t%d\t%.1f\t%.1f\t%.1f\t%.4f\t%.4f\n",
691 i,llen_hist[i],exp(((double)(i))/LN_FACT),
692 score_sums[i],score2_sums[i],mu,var);
693 }
694 }
695 free_hist();
696 fclose(bout);
697 }
698 */
699
s_abort()700 s_abort()
701 {
702 exit(1);
703 }
704