forester/archive/RIO/others/hmmer/src/hmmsearch-pvm.c

   1 /************************************************************
   2  * HMMER - Biological sequence analysis with profile HMMs
   3  * Copyright (C) 1992-1999 Washington University School of Medicine
   4  * All Rights Reserved
   5  *
   6  *     This source code is distributed under the terms of the
   7  *     GNU General Public License. See the files COPYING and LICENSE
   8  *     for details.
   9  ************************************************************/
  10
  11 #ifdef HMMER_PVM
  12
  13 /* hmmsearch-pvm.c
  14  * SRE, Wed Sep 23 09:30:53 1998
  15  *
  16  * PVM slave for hmmsearch.
  17  * RCS $Id: hmmsearch-pvm.c,v 1.1.1.1 2005/03/22 08:34:12 cmzmasek Exp $
  18  */
  19
  20 #include <stdio.h>
  21 #include <stdlib.h>
  22 #include <float.h>
  23 #include <pvm3.h>
  24
  25 #include "version.h"
  26 #include "structs.h"            /* data structures, macros, #define's   */
  27 #include "config.h"             /* compile-time configuration constants */
  28 #include "funcs.h"              /* function declarations                */
  29 #include "globals.h"            /* alphabet global variables            */
  30 #include "squid.h"              /* general sequence analysis library    */
  31
  32 static void leave_pvm(void);
  33
  34 int
  35 main(void)
  36 {
  37   struct plan7_s *hmm;          /* HMM to search with */
  38   struct p7trace_s *tr;         /* trace structure for a Viterbi alignment */
  39   int master_tid;               /* PVM TID of our master */
  40   int alphatype;                /* alphabet type */
  41   int code;                     /* status code for whether we're ok */
  42   int my_idx;                   /* my slave index: 0..nslaves-1, master assigns */
  43   int L;                        /* length of sequence */
  44   char *dsq;                    /* digitized sequence 1..L */
  45   float  sc;                    /* log odds score for seq + HMM */
  46   double pvalue;                /* P-value of sc */
  47   double evalue;                /* bounded E-value of sc (we don't know nseq yet) */
  48   int    do_forward;            /* TRUE to score using Forward() */
  49   int    do_null2;              /* TRUE to use null2 ad hoc correction */
  50   float  globT;                 /* T parameter: keep only hits > globT bits */
  51   double globE;                 /* E parameter: keep hits < globE E-value   */
  52   int    Z;                     /* nseq to base E value calculation on      */
  53   int    nseq;                  /* actual nseq so far (master keeps updating this) */
  54   int    send_trace;            /* TRUE if sc looks significant and we return tr */
  55
  56   /* Register leave_pvm() cleanup function so any exit() call
  57    * first calls pvm_exit().
  58    */
  59   if (atexit(leave_pvm) != 0) { pvm_exit(); Die("slave couldn't register leave_pvm()"); }
  60
  61   /*****************************************************************
  62    * Initialization.
  63    * Master broadcasts the problem to us:
  64    * globT, globE, Z, do_forward, do_null2, alphabet type, HMM,
  65    ******************************************************************/
  66
  67   master_tid = pvm_parent();    /* who's our master? */
  68   my_idx     = -1;
  69
  70   /* wait for a HMMPVM_INIT message, and unpack it;
  71    * get options, set alphabet type, get HMM.
  72    */
  73   pvm_recv(master_tid, HMMPVM_INIT);
  74   pvm_upkfloat(&globT, 1, 1);
  75   pvm_upkdouble(&globE, 1, 1);
  76   pvm_upkint(&Z,          1, 1);
  77   pvm_upkint(&do_forward, 1, 1);
  78   pvm_upkint(&do_null2,   1, 1);
  79   pvm_upkint(&alphatype,  1, 1);
  80   SetAlphabet(alphatype);
  81   hmm = PVMUnpackHMM();
  82
  83   P7Logoddsify(hmm, TRUE);
  84
  85   /* tell the master we're OK and ready to go (or not)
  86    */
  87   code = HMMPVM_OK;
  88   if (hmm == NULL) code = HMMPVM_BAD_INIT;
  89   pvm_initsend(PvmDataDefault);
  90   pvm_pkint(&code, 1, 1);
  91   PVMPackString(RELEASE);
  92   pvm_send(master_tid, HMMPVM_RESULTS);
  93
  94   /*****************************************************************
  95    * Main loop.
  96    * Receive a digitized sequence to search against.
  97    *****************************************************************/
  98
  99   for (;;)
 100     {
 101       SQD_DPRINTF1(("Slave about to do a blocking receive, waiting for input.\n"));
 102       pvm_recv(master_tid, HMMPVM_WORK);
 103       pvm_upkint(&nseq, 1, 1);
 104       if (nseq == -1) break;    /* shutdown signal */
 105       if (my_idx == -1) my_idx = nseq;
 106       pvm_upkint(&L,    1, 1);
 107       SQD_DPRINTF1(("Slave received nseq=%d L=%d my_idx=%d\n", nseq, L, my_idx));
 108       dsq = MallocOrDie(sizeof(char) * (L + 2));
 109       pvm_upkbyte(dsq, L+2, 1);
 110       SQD_DPRINTF1(("Slave unpacked a seq of %d bytes; beginning processing\n", L+2));
 111
 112       /* Score sequence, do alignment (Viterbi), recover trace
 113        */
 114       if (P7ViterbiSize(L, hmm->M) <= RAMLIMIT)
 115         {
 116           SQD_DPRINTF1(("Slave doing Viterbi after estimating %d MB\n", (P7ViterbiSize(L, hmm->M))));
 117           sc = P7Viterbi(dsq, L, hmm, &tr);
 118         }
 119       else
 120         {
 121           SQD_DPRINTF1(("Slave going small after estimating %d MB\n", (P7ViterbiSize(L, hmm->M))));
 122           sc = P7SmallViterbi(dsq, L, hmm, &tr);
 123         }
 124
 125       if (do_forward) sc  = P7Forward(dsq, L, hmm, NULL);
 126       if (do_null2)   sc -= TraceScoreCorrection(hmm, tr, dsq);
 127
 128       pvalue = PValue(hmm, sc);
 129       evalue = Z ? (double) Z * pvalue : (double) nseq * pvalue;
 130       send_trace = (sc >= globT && evalue <= globE) ? 1 : 0;
 131
 132       /* return output
 133        */
 134       SQD_DPRINTF1(("Slave has a result (sc = %.1f); sending back to master\n", sc));
 135       pvm_initsend(PvmDataDefault);
 136       pvm_pkint   (&my_idx,  1, 1);
 137       pvm_pkfloat (&sc,      1, 1);
 138       pvm_pkdouble(&pvalue,  1, 1);
 139       pvm_pkint(&send_trace, 1, 1); /* flag for whether a trace structure is coming */
 140       if (send_trace) PVMPackTrace(tr);
 141       pvm_send(master_tid, HMMPVM_RESULTS);
 142
 143       /* cleanup
 144        */
 145       free(dsq);
 146       P7FreeTrace(tr);
 147     }
 148
 149   /***********************************************
 150    * Cleanup, return.
 151    ***********************************************/
 152
 153   SQD_DPRINTF1(("Slave is done; performing a normal exit.\n"));
 154   FreePlan7(hmm);
 155   exit(0);                      /* pvm_exit() gets called by atexit() registration. */
 156 }
 157
 158 /* Function: leave_pvm()
 159  *
 160  * Purpose:  Cleanup function, to deal with crashes. We register
 161  *           this function using atexit() so it gets called before
 162  *           the slave dies.
 163  */
 164 void leave_pvm(void)
 165 {
 166   SQD_DPRINTF1(("slave leaving PVM.\n"));
 167   pvm_exit();
 168 }
 169
 170
 171 #else /* if HMMER_PVM not defined: include a dummy */
 172
 173 #include <stdio.h>
 174 int main(void)
 175 {
 176   printf("hmmsearch-pvm is disabled. PVM support was not compiled into HMMER.\n");
 177   exit(0);
 178 }
 179
 180 #endif