binaries/src/ViennaRNA/Progs/RNALfold.c

   1 /* Last changed Time-stamp: <2003-04-23 11:56:44 ivo> */
   2 /*
   3                   Ineractive Access to folding Routines
   4
   5                   c Ivo L Hofacker
   6                   Vienna RNA package
   7 */
   8
   9 #include <stdio.h>
  10 #include <stdlib.h>
  11 #include <math.h>
  12 #include <ctype.h>
  13 #include <unistd.h>
  14 #include <string.h>
  15 #include "fold_vars.h"
  16 #include "utils.h"
  17 #include "read_epars.h"
  18 #include "Lfold.h"
  19 #include "RNALfold_cmdl.h"
  20
  21 /*@unused@*/
  22 static char rcsid[] = "$Id: RNALfold.c,v 1.2 2003/07/14 13:38:47 ivo Exp $";
  23
  24 int main(int argc, char *argv[]){
  25   struct  RNALfold_args_info  args_info;
  26   char                        *input_string, *c, *string, *structure, *ParamFile, *ns_bases, *rec_sequence, *rec_id, **rec_rest, *orig_sequence;
  27   int                         i, length, l, sym, r, istty, noconv, maxdist, zsc;
  28   double                      energy, min_en, min_z;
  29   unsigned int                input_type;
  30   unsigned int                rec_type, read_opt;
  31
  32   string        = structure = ParamFile = ns_bases = NULL;
  33   do_backtrack  = 1;
  34   noconv        = 0;
  35   dangles       = 2;
  36   maxdist       = 150;
  37   zsc           = 0;
  38   min_z         = -2.0;
  39   gquad         = 0;
  40   rec_type      = read_opt = 0;
  41   rec_id        = rec_sequence = orig_sequence = NULL;
  42   rec_rest      = NULL;
  43
  44   /*
  45   #############################################
  46   # check the command line parameters
  47   #############################################
  48   */
  49   if(RNALfold_cmdline_parser (argc, argv, &args_info) != 0) exit(1);
  50   /* temperature */
  51   if(args_info.temp_given)        temperature = args_info.temp_arg;
  52   /* do not take special tetra loop energies into account */
  53   if(args_info.noTetra_given)     tetra_loop=0;
  54   /* set dangle model */
  55   if(args_info.dangles_given){
  56     if((args_info.dangles_arg < 0) || (args_info.dangles_arg > 3))
  57       warn_user("required dangle model not implemented, falling back to default dangles=2");
  58     else
  59       dangles = args_info.dangles_arg;
  60   }
  61   /* do not allow weak pairs */
  62   if(args_info.noLP_given)        noLonelyPairs = 1;
  63   /* do not allow wobble pairs (GU) */
  64   if(args_info.noGU_given)        noGU = 1;
  65   /* do not allow weak closing pairs (AU,GU) */
  66   if(args_info.noClosingGU_given) no_closingGU = 1;
  67   /* do not convert DNA nucleotide "T" to appropriate RNA "U" */
  68   if(args_info.noconv_given)      noconv = 1;
  69   /* set energy model */
  70   if(args_info.energyModel_given) energy_set = args_info.energyModel_arg;
  71   /* take another energy parameter set */
  72   if(args_info.paramFile_given)   ParamFile = strdup(args_info.paramFile_arg);
  73   /* Allow other pairs in addition to the usual AU,GC,and GU pairs */
  74   if(args_info.nsp_given)         ns_bases = strdup(args_info.nsp_arg);
  75   /* set the maximum base pair span */
  76   if(args_info.span_given)        maxdist = args_info.span_arg;
  77   if(args_info.zscore_given){
  78 #ifdef USE_SVM
  79     zsc = 1;
  80     if(args_info.zscore_arg != -2)
  81       min_z = args_info.zscore_arg;
  82 #else
  83   nrerror("\'z\' option is available only if compiled with SVM support!");
  84 #endif
  85   }
  86   /* gquadruplex support */
  87   if(args_info.gquad_given)      gquad = 1;
  88
  89   /* check for errorneous parameter options */
  90   if(maxdist < 0){
  91     RNALfold_cmdline_parser_print_help();
  92     exit(EXIT_FAILURE);
  93   }
  94
  95   /* free allocated memory of command line data structure */
  96   RNALfold_cmdline_parser_free (&args_info);
  97
  98   /*
  99   #############################################
 100   # begin initializing
 101   #############################################
 102   */
 103   if (ParamFile != NULL)
 104     read_parameter_file(ParamFile);
 105
 106   if (ns_bases != NULL) {
 107     nonstandards = space(33);
 108     c=ns_bases;
 109     i=sym=0;
 110     if (*c=='-') {
 111       sym=1; c++;
 112     }
 113     while (*c!='\0') {
 114       if (*c!=',') {
 115         nonstandards[i++]=*c++;
 116         nonstandards[i++]=*c;
 117         if ((sym)&&(*c!=*(c-1))) {
 118           nonstandards[i++]=*c;
 119           nonstandards[i++]=*(c-1);
 120         }
 121       }
 122       c++;
 123     }
 124   }
 125
 126   istty = isatty(fileno(stdout))&&isatty(fileno(stdin));
 127   read_opt |= VRNA_INPUT_NO_REST;
 128   if(istty){
 129     print_tty_input_seq();
 130     read_opt |= VRNA_INPUT_NOSKIP_BLANK_LINES;
 131   }
 132   /*
 133   #############################################
 134   # main loop: continue until end of file
 135   #############################################
 136   */
 137   while(
 138     !((rec_type = read_record(&rec_id, &rec_sequence, &rec_rest, read_opt))
 139         & (VRNA_INPUT_ERROR | VRNA_INPUT_QUIT))){
 140
 141     /*
 142     ########################################################
 143     # init everything according to the data we've read
 144     ########################################################
 145     */
 146     if(rec_id && !istty) printf("%s\n", rec_id);
 147
 148     length = (int)strlen(rec_sequence);
 149
 150     /* convert DNA alphabet to RNA if not explicitely switched off */
 151     if(!noconv) str_DNA2RNA(rec_sequence);
 152     /* store case-unmodified sequence */
 153     orig_sequence = strdup(rec_sequence);
 154     /* convert sequence to uppercase letters only */
 155     str_uppercase(rec_sequence);
 156
 157     if(istty) printf("length = %d\n", length);
 158     /*
 159     ########################################################
 160     # done with 'stdin' handling
 161     ########################################################
 162     */
 163
 164     min_en = (zsc) ? Lfoldz((const char *)rec_sequence, NULL, maxdist, zsc, min_z) : Lfold((const char *)rec_sequence, NULL, maxdist);
 165     printf("%s\n", orig_sequence);
 166
 167     if (istty)
 168       printf("\n minimum free energy = %6.2f kcal/mol\n", min_en);
 169     else
 170       printf(" (%6.2f)\n", min_en);
 171
 172     (void) fflush(stdout);
 173
 174     /* clean up */
 175     if(rec_id) free(rec_id);
 176     free(rec_sequence);
 177     free(orig_sequence);
 178     rec_id = rec_sequence = orig_sequence = NULL;
 179     rec_rest = NULL;
 180     /* print user help for the next round if we get input from tty */
 181
 182     if(istty) print_tty_input_seq();
 183   }
 184   return EXIT_SUCCESS;
 185 }
 186