2 /* copyright (c) 1996, 1997, 1998, 1999 William R. Pearson and the
5 /* $Name: fa_34_26_5 $ - $Id: mshowbest.c,v 1.44 2006/06/30 19:46:36 wrp Exp $ */
7 /* 29-Oct-2003 - changes so that bbp->cont < 0 => aa1 sequence is
8 already in aa1, no re_openlib or re_getlib required
11 /* 14-May-2003 Changes to use a more consistent coordinate numbering
12 system for displays. aln->d_start[01] is now consistently used
13 to report the start of the alignment in all functions, and
14 mshowbest.c has been modified to use d_start[01] instead of
15 d_start[01]-1. aln->min[01] now starts at 0 for all functions;
16 instead of 1 for some functions (dropnfa.c, dropgsw.c, dropfs2.c
39 /* function calls necessary to re_getlib() the sequence and, do
40 alignments, if necessary
43 #define RANLIB (m_fptr->ranlib)
46 re_getlib(unsigned char *, int, int, int, int, int, long *, long *,
47 struct lmf_str *m_fptr);
49 #include "drop_func.h"
51 struct lmf_str *re_openlib(struct lmf_str *, int outtty);
54 extern void cal_coord(int n0, int n1, long sq0off, long loffset,
55 struct a_struct *aln);
57 void header_aux(FILE *);
58 void show_aux(FILE *, struct beststr *);
59 void w_abort (char *p, char *p1);
61 /* BBP_INFO get stuff directly from beststr or from beststr->desptr */
63 #define BBP_INFO(info) bbp->desptr->info
65 #define BBP_INFO(info) bbp->info
68 extern double zs_to_bit(double, int, int);
70 /* showbest() shows a list of high scoring sequence descriptions, and
71 their scores. If -m 9, then an additional complete set of
72 alignment information is provided.
74 If PCOMPLIB or m_msg.quiet then the number of high scores to be
75 shown is pre-determined by m_msg.mshow before showbest is called.
77 The comp_lib.c version re_getlib()'s the sequence for its
78 discription, and then does another alignment for -m 9 (Thus, it
79 needs an f_str. The PCOMPLIB version has everything available in
80 beststr before showbest() is called.
83 void showbest (FILE *fp,
85 unsigned char **aa0, unsigned char *aa1, int maxn,
87 struct beststr **bptr,int nbest, int qlib, struct mngmsg *m_msg,
88 struct pstruct pst, struct db_str db,
96 char bline[MAX_BLINE], fmt[40], pad[MAX_BLINE], rline[40];
98 int istart = 0, istop, ib;
106 char tmp_str[20], *seqc;
111 int lc, maxc, nident, ngap;
112 float percent, gpercent;
113 struct a_struct *aln_p;
118 struct lmf_str *m_fptr;
121 strncpy(rel_label,"\0",2);
123 strncpy(rel_label," related",sizeof(rel_label));
126 strncpy(rel_label," unrelated",sizeof(rel_label));
128 rel_label[sizeof(rel_label)-1]='\0';
133 quiet = m_msg->quiet;
138 if (m_msg->aln.llen > MAX_BLINE) m_msg->aln.llen = MAX_BLINE;
140 if (pst.zsflag < 0) r_margin = 10;
141 else if (pst.zsflag>=0 && m_msg->srelv > 1 ) r_margin = 19;
144 if (m_msg->markx & MX_M9SUMM && m_msg->show_code == SHOW_CODE_ID) {
152 if (m_msg->nframe < 0)
154 sprintf(fmt,"%%-%ds (%%4d)",m_msg->aln.llen-r_margin);
156 sprintf(fmt,"%%-%ds [%%4d](%%4d)",m_msg->aln.llen-(r_margin+4));
159 sprintf(fmt,"%%-%ds (%%4d)",m_msg->aln.llen-(r_margin+4));
161 memset(pad,' ',m_msg->aln.llen-(r_margin+6));
162 pad[m_msg->aln.llen-(r_margin+12)]='\0';
164 if (quiet != -1) { /* quiet is set to -1 in comp_mlib.c to force
165 all significant hits to be shown */
167 if (m_msg->mshow == -1) nshow = nbest; /* show all */
168 /* show specified number */
169 else if (m_msg->mshow_flg) {
170 nshow = min (m_msg->mshow, nshow);
173 else nshow = m_msg->nshow;
175 if (quiet==0) istop = 20;
179 printf(" How many scores would you like to see? [%d] ",m_msg->nshow);
181 if (fgets(rline,20,stdin)==NULL) exit(0);
182 nshow = m_msg->nshow;
183 if (rline[0]!='\n' && rline[0]!=0) sscanf(rline,"%d",&nshow);
184 if (nshow<=0) nshow = min(20,nbest);
187 if ((bp = strchr (m_msg->qtitle, '\n')) != NULL) *bp = '\0';
188 /* fprintf (fp, "%3d %s\n", qlib,m_msg->qtitle); */
190 if (m_msg->markx & MX_HTML) fprintf(fp,"<p><tt><pre>\n");
192 if (pst.zsflag >= 0) {
193 if (bptr[0]->escore < m_msg->e_cut) {
194 if (m_msg->z_bits==1) {/* show bit score */
195 fprintf(fp,"The best%s scores are:%s%s bits E(%ld)",
196 rel_label,pad,m_msg->label,pst.zdb_size);
198 else {/* show z-score */
199 fprintf(fp,"The best%s scores are:%s%s z-sc E(%ld)",
200 rel_label,pad,m_msg->label,pst.zdb_size);
203 if (m_msg->markx & MX_M9SUMM) {
204 if (m_msg->show_code == SHOW_CODE_ID) {
206 fprintf(fp," %%_id %%_sim alen");
208 fprintf(fp," %%_id alen");
212 if (m_msg->markx & MX_HTML && m_msg->show_code !=1) { fprintf(fp,"<!-- ");}
214 fprintf(fp,"\t%%_id %%_gid %4s alen an0 ax0 pn0 px0 an1 ax1 pn1 px1 gapq gapl fs ",m_msg->f_id1);
216 fprintf(fp,"\t%%_id %%_sim %4s alen an0 ax0 pn0 px0 an1 ax1 pn1 px1 gapq gapl fs ",m_msg->f_id1);
219 if (m_msg->show_code == SHOW_CODE_ALIGN) { fprintf(fp," aln_code"); }
220 if (m_msg->markx & MX_HTML && m_msg->show_code!=1) { fprintf(fp," -->");}
225 fprintf(fp,"!! No library sequences with E() < %.2g\n",m_msg->e_cut);
227 if (m_msg->markx & MX_HTML) fprintf(fp,"<p></tt></pre>\n");
232 fprintf(fp,"The best%s scores are:%s%s",rel_label,pad,m_msg->label);
234 if (m_msg->markx & MX_M9SUMM) {
235 if (m_msg->show_code == SHOW_CODE_ID) {
237 fprintf(fp," %%_id %%_sm alen");
239 fprintf(fp," %%_id alen");
244 fprintf(fp,"\t%%_id %%_gid %4s alen an0 ax0 pn0 px0 an1 ax1 pn1 px1 gapq gapl fs ",m_msg->f_id1);
246 fprintf(fp,"\t%%_id %%_sim %4s alen an0 ax0 pn0 px0 an1 ax1 pn1 px1 gapq gapl fs ",m_msg->f_id1);
250 if (m_msg->show_code == SHOW_CODE_ALIGN) { fprintf(fp," aln_code"); }
256 istop = min(nbest,nshow);
257 for (ib=istart; ib<istop; ib++) {
260 if (BBP_INFO(nsfnum) > 0 && sfn_cmp(m_msg->qsfnum_n,BBP_INFO(sfnum))) continue;
262 if (BBP_INFO(nsfnum) > 0 && sfn_cmp(m_msg->qsfnum,BBP_INFO(sfnum))) {
263 istop = min(istop+1,nbest);
265 fprintf(stderr,"skipping %d: %d==%d\n",ib,m_msg->qsfnum,BBP_INFO(sfnum));
271 if (BBP_INFO(nsfnum) == 0 || (BBP_INFO(nsfnum) > 0 && !sfn_cmp(m_msg->qsfnum,BBP_INFO(sfnum)))) {
272 istop = min(istop+1,nbest);
274 fprintf(stderr,"skipping %d: %d==%d\n",ib,m_msg->qsfnum,BBP_INFO(sfnum));
280 if (quiet==1 && pst.zsflag>=0) {
281 if (bbp->escore > m_msg->e_cut) {
285 else if (bbp->escore < m_msg->e_low) continue;
289 if ((m_fptr=re_openlib(bbp->m_file_p,!m_msg->quiet))==NULL) {
290 fprintf(stderr,"*** cannot re-open %s\n",bbp->m_file_p->lb_name);
293 RANLIB(bline,m_msg->aln.llen,bbp->lseek,bbp->libstr,m_fptr);
295 strncpy(bline,BBP_INFO(bline),m_msg->aln.llen-r_margin);
296 bline[m_msg->aln.llen]='\0';
299 /* l_name is used to build an HTML link from the bestscore line to
300 the alignment. It can also be used to discriminate multiple hits
301 from the same long sequence. This requires that fast_pan use -m 6. */
303 strncpy(l_name,bline,sizeof(l_name)); /* get rid of text after second "|" */
304 l_name[sizeof(l_name)-1]='\0';
305 if ((bp=strchr(l_name,' '))!=NULL) *bp=0;
306 if ((bp=strchr(&l_name[3],'|'))!=NULL) *bp='\0';
307 if (m_msg->nframe > 2) sprintf(&l_name[strlen(l_name)],"_%d",bbp->frame+1);
308 else if (m_msg->nframe > 0 && bbp->frame == 1)
309 strncat(l_name,"_r",sizeof(l_name));
310 if (bbp->cont-1 > 0) {
311 sprintf(tmp_str,":%d",bbp->cont-1);
312 strncat(l_name,tmp_str,sizeof(l_name)-strlen(l_name));
317 if (m_msg->stages>1 || m_msg->markx & MX_M9SUMM) {
318 if (bbp->cont >= 0) {
319 n1 = re_getlib(aa1,maxn,m_msg->maxt3,m_msg->loff,bbp->cont,m_msg->term_code,
320 &loffset,&l_off,bbp->m_file_p);
323 if (! m_msg->markx & MX_M9SUMM) {
324 do_opt (aa0[bbp->frame], m_msg->n0, aa1, n1, bbp->frame, &pst, f_str[bbp->frame], &rst);
325 bbp->score[2]=rst.score[2];
329 do_walign(aa0[bbp->frame],m_msg->n0, aa1, n1, bbp->frame,
330 &pst, f_str[bbp->frame], &bbp->a_res, &bbp->have_ares);
333 /* save the alignment encoding for future use */
334 if (bbp->have_ares && ((tres = calloc(bbp->a_res.nres+1,sizeof(int)))!=NULL)) {
335 memcpy(tres,bbp->a_res.res,sizeof(int)*bbp->a_res.nres);
336 bbp->a_res.res = tres;
339 aln_func_vals(bbp->frame, &m_msg->aln);
341 maxc = bbp->a_res.nres + 4*m_msg->aln.llen+4;
344 if (m_msg->show_code == SHOW_CODE_ALIGN) {
345 if ((seqc=(char *)calloc(maxc,sizeof(char)))!=NULL) {
346 lc=calc_code(aa0[bbp->frame],m_msg->n0,
348 &m_msg->aln,bbp->a_res,
349 pst,seqc,maxc,f_str[bbp->frame]);
350 seqc_len = strlen(seqc);
354 lc=calc_id(aa0[bbp->frame],m_msg->n0,aa1,n1,
355 &m_msg->aln, bbp->a_res,
356 pst,f_str[bbp->frame]);
358 m_msg->aln.a_len = lc;
360 nident = m_msg->aln.nident;
361 if (lc > 0) percent = (100.0*(float)nident)/(float)lc;
362 else percent = -1.00;
364 ngap = m_msg->aln.ngap_q + m_msg->aln.ngap_l;
366 if (lc-ngap > 0) gpercent = (100.0*(float)nident)/(float)(lc-ngap);
367 else gpercent = -1.00;
369 if (lc-ngap > 0) gpercent = (100.0*(float)m_msg->aln.nsim)/(float)(lc);
370 else gpercent = -1.00;
377 n1tot = (BBP_INFO(n1tot_p)) ? *BBP_INFO(n1tot_p) : bbp->n1;
380 if ((m_msg->markx & MX_HTML) && !strncmp(bline,"gi|",3)) {
381 bp = strchr(bline+4,'|')+1;
383 gi_num = atoi(bline+3);
387 bp[m_msg->aln.llen-r_margin]='\0';
389 bp[m_msg->aln.llen-r_margin-5]='\0';
392 if (m_msg->nframe == -1) bp[m_msg->aln.llen-r_margin]='\0';
393 else bp[m_msg->aln.llen-(r_margin+4)]='\0';
396 fprintf (fp, fmt,bp,n1tot);
398 if (m_msg->nframe == -1) {
399 fprintf (fp, fmt,bp,BBP_INFO(sfnum[0]),n1tot);
401 else {fprintf (fp, fmt,bp,n1tot);}
404 if (m_msg->nframe > 2) fprintf (fp, " [%d]", bbp->frame+1);
405 else if (m_msg->nframe >= 0) fprintf(fp," [%c]",(bbp->frame > 0 ?'r':'f'));
407 if (m_msg->srelv == 1) fprintf (fp, " %4d", bbp->score[pst.score_ix]);
409 if (m_msg->srelv-1 > 0) fprintf (fp, " %4d", bbp->score[0]);
410 if (m_msg->srelv-1 > 1 || m_msg->stages>1)
411 fprintf (fp, " %4d", bbp->score[1]);
412 fprintf (fp, " %4d", bbp->score[pst.score_ix]);
416 if (m_msg->z_bits==1) {
417 fprintf (fp, " %.1f %7.2g",zs_to_bit(bbp->zscore,m_msg->n0,bbp->n1),bbp->escore);
419 else fprintf (fp, " %.1f %7.2g",bbp->zscore,bbp->escore);
425 percent = bbp->percent;
426 gpercent = bbp->gpercent;
428 seqc = bbp->aln_code;
429 seqc_len = bbp->aln_code_n;
430 loffset = bbp->desptr->loffset;
433 aln_p = &(m_msg->aln);
436 if (m_msg->markx & MX_M9SUMM) {
437 if (m_msg->show_code != SHOW_CODE_ID) {
438 if (m_msg->markx & MX_HTML) fprintf(fp,"<!-- ");
439 cal_coord(m_msg->n0,bbp->n1,m_msg->sq0off,loffset+l_off-1,aln_p);
441 /* %_id %_sim s-w alen an0 ax0 pn0 px0 an1 ax1 pn1 px1 gapq gapl fs */
442 /* alignment min max min max */
443 /* sequence coordinate min max min max */
444 fprintf(fp,"\t%5.3f %5.3f %4d %4d %4ld %4ld %4ld %4ld %4ld %4ld %4ld %4ld %3d %3d %3d",
445 percent/100.0,gpercent/100.0, bbp->sw_score,aln_p->a_len,
446 aln_p->d_start0,aln_p->d_stop0,
447 m_msg->sq0off, m_msg->sq0off+m_msg->n0-1,
448 aln_p->d_start1,aln_p->d_stop1,
449 loffset+l_off, loffset+l_off+bbp->n1-1,
450 aln_p->ngap_q,aln_p->ngap_l,aln_p->nfs);
451 if (m_msg->show_code == SHOW_CODE_ALIGN
452 && seqc_len > 0 && seqc != NULL) {
453 fprintf(fp,"\t%s",seqc);
454 /* fprintf(fp," [%2d:%d]",bbp->wrkr,bbp->seqnm); */
458 if (m_msg->markx & MX_HTML) fprintf(fp," -->");
462 fprintf(fp," %5.3f %5.3f %4d", percent/100.0,(float)aln_p->nsim/(float)aln_p->a_len,aln_p->a_len);
464 fprintf(fp," %5.3f %4d", percent/100.0,aln_p->a_len);
468 if (m_msg->markx & MX_HTML) fprintf(fp," <A HREF=\"#%s\">align</A>",l_name);
474 printf(" More scores? [0] ");
476 if (fgets(rline,20,stdin)==NULL) exit(0);
478 if (rline[0]!='\n' && rline[0]!=0) sscanf(rline,"%d",&ntmp);
479 if (ntmp<=0) ntmp = 0;
487 if (ib < nbest && (pst.zsflag>=0 && bbp->escore < m_msg->e_cut)) {
488 if (m_msg->mshow_flg && istop >= m_msg->mshow) goto done;
495 m_msg->nshow = nshow;
496 if (m_msg->markx & MX_HTML) fprintf(fp,"</pre></tt><p><hr><p>\n");
497 if (fp!=stdout) fprintf(fp,"\n");
501 q[] has one set of sfnums, 0 terminated
503 return first match or 0