mshowbest.c

序列对齐 Compare a protein sequence to a protein sequence database or a DNA sequence to a DNA sequenc

/* copyright (c) 1996, 1997, 1998, 1999 William R. Pearson and the
   U. of Virginia */

/* $Name: fa35_03_06 $ - $Id: mshowbest.c,v 1.54 2008/01/11 15:27:43 wrp Exp $ */

/*   29-Oct-2003 - changes so that bbp->seq->cont < 0 => aa1 sequence is
     already in aa1, no re_openlib or re_getlib required
*/

/*   14-May-2003 Changes to use a more consistent coordinate numbering
     system for displays.  aln->d_start[01] is now consistently used
     to report the start of the alignment in all functions, and
     mshowbest.c has been modified to use d_start[01] instead of
     d_start[01]-1.  aln->min[01] now starts at 0 for all functions;
     instead of 1 for some functions (dropnfa.c, dropgsw.c, dropfs2.c
     earlier).
*/

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#include "defs.h"
#include "structs.h"
#include "param.h"

#ifndef PCOMPLIB
#include "mm_file.h"
#include "mw.h"
#else
#include "p_mw.h"
#endif


#define MAX_BLINE 256

#ifndef PCOMPLIB
/* function calls necessary to re_getlib() the sequence and, do
   alignments, if necessary
*/

#define RANLIB (m_fptr->ranlib)

int
re_getlib(unsigned char *, int, int, int, int, int, long *, long *, 
	  struct lmf_str *m_fptr);

#include "drop_func.h"

struct lmf_str *re_openlib(struct lmf_str *, int outtty);
#endif

extern void cal_coord(int n0, int n1, long qoffset, long loffset,
		      struct a_struct *aln);

void header_aux(FILE *);
void show_aux(FILE *, struct beststr *);
void w_abort (char *p, char *p1);

/* BBP_INFO get stuff directly from beststr or from beststr->desptr */
#ifdef PCOMPLIB
#define BBP_INFO(info) bbp->seq->info
#else
#define BBP_INFO(info) bbp->seq->info
#endif

extern double zs_to_bit(double, int, int);

/* showbest() shows a list of high scoring sequence descriptions, and
   their rst.scores.  If -m 9, then an additional complete set of
   alignment information is provided.

   If PCOMPLIB or m_msg.quiet then the number of high scores to be
   shown is pre-determined by m_msg.mshow before showbest is called.

   The comp_lib2.c version re_getlib()'s the sequence for its
   discription, and then does another alignment for -m 9 (Thus, it
   needs an f_str.  The PCOMPLIB version has everything available in
   beststr before showbest() is called.
*/

void showbest (FILE *fp, 
#ifndef PCOMPLIB
	       unsigned char **aa0, unsigned char *aa1save, int maxn,
#endif
	       struct beststr **bptr,int nbest,
	       int qlib, struct mngmsg *m_msg,
	       struct pstruct *ppst, struct db_str db,
	       char **info_gstring2
#ifndef PCOMPLIB
	       ,void **f_str
#endif
)
{
  unsigned char *aa1;
  int ntmp = 0;
  char bline[MAX_BLINE], fmt[40], pad[MAX_BLINE], rline[40];
  char l_name[128];
  int istart = 0, istop, ib;
  int nshow;
  int quiet;
  int r_margin;
  struct beststr *bbp;
  int n1tot;
  char *bp;
  char rel_label[12];
  char tmp_str[20], *seq_code, *ann_code;
  int seq_code_len, ann_code_len;
  long loffset;		/* loffset is offset from beginning of real sequence */
  long l_off;		/* l_off is the the virtual coordinate of residue 1 */
  int n0, n1;
  struct rstruct rst;
  int lc, seqc_max, annc_max, nident, ngap;
  float percent, gpercent;
  struct a_struct *aln_p;
  int *tres;
  int gi_num;

#ifndef PCOMPLIB
  struct lmf_str *m_fptr;
#endif

  strncpy(rel_label,"\0",2);
#ifdef SHOWREL
  strncpy(rel_label," related",sizeof(rel_label));
#endif
#ifdef SHOWUN
  strncpy(rel_label," unrelated",sizeof(rel_label));
#endif
  rel_label[sizeof(rel_label)-1]='\0';

#ifdef PCOMPLIB
  quiet = 1;
#else
  quiet = m_msg->quiet;
#endif

  n0 = m_msg->n0;

  if (m_msg->aln.llen > MAX_BLINE) m_msg->aln.llen = MAX_BLINE;

  if (ppst->zsflag < 0) r_margin = 10;
  else if (ppst->zsflag>=0  && m_msg->srelv > 1 ) r_margin = 19;
  else r_margin = 10;

  if (m_msg->markx & MX_M9SUMM && m_msg->show_code == SHOW_CODE_ID) {
#ifdef SHOWSIM
    r_margin += 15;
#else
    r_margin += 10;
#endif
  }

  if (m_msg->nframe < 0) {
#ifndef SUPERFAMNUM
    sprintf(fmt,"%%-%ds (%%4d)",m_msg->aln.llen-r_margin);
#else
    sprintf(fmt,"%%-%ds [%%4d](%%4d)",m_msg->aln.llen-(r_margin+4));
#endif
  }
  else { sprintf(fmt,"%%-%ds (%%4d)",m_msg->aln.llen-(r_margin+4)); }

  memset(pad,' ',m_msg->aln.llen-(r_margin+6));
  pad[m_msg->aln.llen-(r_margin+12)]='\0';

  if (quiet != -1) {	/* quiet is set to -1 in comp_lib2.c to force
			   all significant hits to be shown */
    nshow = 20;
    if (m_msg->mshow == -1) {nshow = nbest;}		/* show all */
    /* show specified number */
    else if (m_msg->mshow_flg) {
      nshow = min (m_msg->mshow, nshow);
    }
  }
  else nshow = m_msg->nshow;

  if (quiet==0) istop = 20;
  else istop = nshow;

  if (quiet==0) {
    printf(" How many scores would you like to see? [%d] ",m_msg->nshow);
    fflush(stdout);
    if (fgets(rline,20,stdin)==NULL) exit(0);
    nshow = m_msg->nshow;
    if (rline[0]!='\n' && rline[0]!=0) sscanf(rline,"%d",&nshow);
    if (nshow<=0) nshow = min(20,nbest);
  }

  if ((bp = strchr (m_msg->qtitle, '\n')) != NULL) *bp = '\0';
/*   fprintf (fp, "%3d %s\n", qlib,m_msg->qtitle); */

  if (m_msg->markx & MX_HTML) fprintf(fp,"<p><tt><pre>\n");

  if (ppst->zsflag >= 0) {
    if (bptr[0]->rst.escore < m_msg->e_cut) {
      if (m_msg->z_bits==1) {/* show bit score */
	fprintf(fp,"\nThe best%s scores are:%s%s bits E(%ld)",
	       rel_label,pad,m_msg->label,ppst->zdb_size);
      }
      else {/* show z-score */
	fprintf(fp,"\nThe best%s scores are:%s%s z-sc E(%ld)",
	      rel_label,pad,m_msg->label,ppst->zdb_size);
      }
      header_aux(fp);
      if (m_msg->markx & MX_M9SUMM) {
	if (m_msg->show_code == SHOW_CODE_ID) {
#ifdef SHOWSIM
	  fprintf(fp," %%_id  %%_sim  alen");
#else
	  fprintf(fp," %%_id  alen");
#endif
	}
	else {
	if (m_msg->markx & MX_HTML && m_msg->show_code !=1) { fprintf(fp,"<!-- ");}
#ifndef SHOWSIM
	  fprintf(fp,"\t%%_id  %%_gid %4s  alen  an0  ax0  pn0  px0  an1  ax1 pn1 px1 gapq gapl  fs ",m_msg->f_id1);
#else
	  fprintf(fp,"\t%%_id  %%_sim %4s  alen  an0  ax0  pn0  px0  an1  ax1 pn1 px1 gapq gapl  fs ",m_msg->f_id1);
#endif
	}
	if (m_msg->show_code == SHOW_CODE_ALIGN) {	fprintf(fp," aln_code"); }
	if (m_msg->markx & MX_HTML && m_msg->show_code!=1) { fprintf(fp," -->");}
      }
      fprintf(fp,"\n");
    }
    else {
      fprintf(fp,"!! No library sequences with E() < %.2g\n",m_msg->e_cut);
      m_msg->nshow = 0;
      if (m_msg->markx & MX_HTML) fprintf(fp,"<p></tt></pre>\n");
      return;
    }
  }
  else {
    fprintf(fp,"\nThe best%s scores are:%s%s",rel_label,pad,m_msg->label);
    header_aux(fp);
    if (m_msg->markx & MX_M9SUMM) {
      if (m_msg->show_code == SHOW_CODE_ID) {
#ifdef SHOWSIM
	fprintf(fp," %%_id  %%_sm  alen");
#else
	fprintf(fp," %%_id  alen");
#endif
      }
      else {
#ifndef SHOWSIM
	fprintf(fp,"\t%%_id  %%_gid %4s  alen  an0  ax0  pn0  px0  an1  ax1 pn1 px1 gapq gapl  fs ",m_msg->f_id1);
#else
	fprintf(fp,"\t%%_id  %%_sim %4s  alen  an0  ax0  pn0  px0  an1  ax1 pn1 px1 gapq gapl  fs ",m_msg->f_id1);
#endif	/* SHOWSIM */
      }
    }
    if (m_msg->show_code == SHOW_CODE_ALIGN) {	fprintf(fp," aln_code"); }
    fprintf(fp,"\n");
  }

  istart = 0;
l1:
  istop = min(nbest,nshow);
  for (ib=istart; ib<istop; ib++) {
    bbp = bptr[ib];

#ifndef PCOMPLIB
#ifdef DEBUG
    if (bbp->seq->n1 != bbp->n1 ) {
      fprintf(stderr, " *** lib len error [%d!=%d] *** %s score %d\n",
	      bbp->seq->n1,bbp->n1, bbp->seq->libstr, bbp->rst.score[0]);
    }
#endif
#endif

#ifdef SUPERFAMNUM
    if (BBP_INFO(nsfnum) > 0 && sfn_cmp(m_msg->qsfnum_n,BBP_INFO(sfnum))) continue;
#ifdef SHOWUN
    if (BBP_INFO(nsfnum) > 0 && sfn_cmp(m_msg->qsfnum,BBP_INFO(sfnum))) {
      istop = min(istop+1,nbest);
    /*
      fprintf(stderr,"skipping %d: %d==%d\n",ib,m_msg->qsfnum,BBP_INFO(sfnum));
      */
      continue;
    }
#endif	/* SHOWUN */
#ifdef SHOWREL
    if (BBP_INFO(nsfnum) == 0 || (BBP_INFO(nsfnum) > 0 && !sfn_cmp(m_msg->qsfnum,BBP_INFO(sfnum)))) {
      istop = min(istop+1,nbest);
      continue;
    }
#endif	/* SHOWREL */
#endif	/* SUPERFAMNUM */
    if (quiet==1 && ppst->zsflag>=0) {
      if (bbp->rst.escore > m_msg->e_cut) {
	nshow = ib;
	goto done;
      }
      else if (bbp->rst.escore < m_msg->e_low) continue;
    }

#ifndef PCOMPLIB
    if ((m_fptr=re_openlib(bbp->seq->m_file_p,!m_msg->quiet))==NULL) {
      fprintf(stderr,"*** cannot re-open %s\n",bbp->seq->m_file_p->lb_name);
      exit(1);
    }
    RANLIB(bline,m_msg->aln.llen,bbp->seq->lseek,bbp->seq->libstr,m_fptr);
#else	/* PCOMPLIB */
  strncpy(bline,BBP_INFO(bline),m_msg->aln.llen-r_margin);
  bline[m_msg->aln.llen]='\0';
#endif

  /* l_name is used to build an HTML link from the bestscore line to
     the alignment.  It can also be used to discriminate multiple hits
     from the same long sequence.  This requires that fast_pan use -m 6. */