WSTester updated to work plus hopefully all the other changes that need to go into...

[jabaws.git] / binaries / src / ViennaRNA / Progs / RNAaliduplex.c

diff --git a/binaries/src/ViennaRNA/Progs/RNAaliduplex.c b/binaries/src/ViennaRNA/Progs/RNAaliduplex.c
new file mode 100644 (file)
index 0000000..ed0afb2
--- /dev/null
+++ b/binaries/src/ViennaRNA/Progs/RNAaliduplex.c
@@ -0,0 +1,172 @@
+/* Last changed Time-stamp: <2005-07-23 16:50:24 ivo> */
+/*
+             Compute duplex structure of two RNA strands
+
+                           c Ivo L Hofacker
+                          Vienna RNA package
+*/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <ctype.h>
+#include <unistd.h>
+#include <string.h>
+#include "fold.h"
+#include "fold_vars.h"
+#include "utils.h"
+#include "aln_util.h"
+#include "read_epars.h"
+#include "subopt.h"
+#include "duplex.h"
+#include "RNAaliduplex_cmdl.h"
+
+/*@unused@*/
+static char rcsid[] = "$Id: RNAaliduplex.c,v 1.1 2007/08/26 10:08:44 ivo Exp $";
+
+PRIVATE void  print_struc(duplexT const *dup);
+
+#define MAX_NUM_NAMES    500
+
+int main(int argc, char *argv[]){
+  struct        RNAaliduplex_args_info args_info;
+  char          *AS1[MAX_NUM_NAMES], *AS2[MAX_NUM_NAMES], *names1[MAX_NUM_NAMES], *names2[MAX_NUM_NAMES];
+  char          *ParamFile=NULL, *c, *ns_bases=NULL;
+  int           i, r, n_seq, n_seq2;
+  int           istty, delta=-1, sym;
+  int           noconv=0;
+  duplexT       mfe, *subopt;
+  FILE          *file1=NULL, *file2=NULL; /* input alignments */
+
+  dangles = 2;
+
+  /*
+  #############################################
+  # check the command line parameters
+  #############################################
+  */
+  if(RNAaliduplex_cmdline_parser (argc, argv, &args_info) != 0) exit(1);
+  /* temperature */
+  if(args_info.temp_given)        temperature = args_info.temp_arg;
+  /* do not take special tetra loop energies into account */
+  if(args_info.noTetra_given)     tetra_loop=0;
+  /* set dangle model */
+  if(args_info.dangles_given){
+    if((args_info.dangles_arg < 0) || (args_info.dangles_arg > 3))
+      warn_user("required dangle model not implemented, falling back to default dangles=2");
+    else
+     dangles = args_info.dangles_arg;
+  }
+  /* do not allow weak pairs */
+  if(args_info.noLP_given)        noLonelyPairs = 1;
+  /* do not allow wobble pairs (GU) */
+  if(args_info.noGU_given)        noGU = 1;
+  /* do not allow weak closing pairs (AU,GU) */
+  if(args_info.noClosingGU_given) no_closingGU = 1;
+  /* do not convert DNA nucleotide "T" to appropriate RNA "U" */
+  if(args_info.noconv_given)      noconv = 1;
+  /* take another energy parameter set */
+  if(args_info.paramFile_given)   ParamFile = strdup(args_info.paramFile_arg);
+  /* Allow other pairs in addition to the usual AU,GC,and GU pairs */
+  if(args_info.nsp_given)         ns_bases = strdup(args_info.nsp_arg);
+  /*energy range */
+  if(args_info.deltaEnergy_given) delta = (int) (0.1+args_info.deltaEnergy_arg*100);
+  /* sorted output */
+  if(args_info.sorted_given)      subopt_sorted = 1;
+
+  /* check unnamed options a.k.a. filenames of input alignments */
+  if(args_info.inputs_num == 2){
+    file1 = fopen(args_info.inputs[0], "r");
+    if(file1 == NULL){
+      fprintf(stderr, "can't open %s\n", args_info.inputs[0]);
+    }
+    file2 = fopen(args_info.inputs[1], "r");
+    if(file2 == NULL) {
+      fprintf(stderr, "can't open %s\n", args_info.inputs[1]);
+    }
+  }
+  else{
+    RNAaliduplex_cmdline_parser_print_help();
+    exit(1);
+  }
+
+  if (!(file1 && file2)){
+    fprintf(stderr, "No input files");
+    RNAaliduplex_cmdline_parser_print_help();
+    exit(1);
+  }
+  else{
+    n_seq   = read_clustal(file1, AS1, names1);
+    n_seq2  = read_clustal(file2, AS2, names2);
+    fclose(file1);
+    fclose(file2);
+    if(n_seq != n_seq2) nrerror("unequal number of seqs in alignments");
+  }
+
+  /* free allocated memory of command line data structure */
+  RNAaliduplex_cmdline_parser_free (&args_info);
+
+  /*
+  #############################################
+  # begin initializing
+  #############################################
+  */
+
+
+  if (ParamFile != NULL)
+    read_parameter_file(ParamFile);
+
+  if (ns_bases != NULL) {
+    nonstandards = space(33);
+    c=ns_bases;
+    i=sym=0;
+    if (*c=='-') {
+      sym=1; c++;
+    }
+    while (*c!='\0') {
+      if (*c!=',') {
+        nonstandards[i++]=*c++;
+        nonstandards[i++]=*c;
+        if ((sym)&&(*c!=*(c-1))) {
+          nonstandards[i++]=*c;
+          nonstandards[i++]=*(c-1);
+        }
+      }
+      c++;
+    }
+  }
+
+  istty = isatty(fileno(stdout))&&isatty(fileno(stdin));
+
+  /*
+  #############################################
+  # begin calculations
+  #############################################
+  */
+  update_fold_params();
+
+  if (delta>=0) {
+    duplexT *sub;
+    subopt = aliduplex_subopt((const char **)AS1, (const char **)AS2, delta, 5);
+    for (sub=subopt; sub->i >0; sub++) {
+      print_struc(sub);
+      free(sub->structure);
+    }
+    free(subopt);
+  } else {
+    mfe = aliduplexfold((const char **)AS1, (const char **)AS2);
+    print_struc(&mfe);
+    free(mfe.structure);
+  }
+  for (i=0; AS1[i]; i++) {
+    free(AS1[i]); free(names1[i]);
+    free(AS2[i]); free(names2[i]);
+  }
+  return 0;
+}
+
+PRIVATE void print_struc(duplexT const *dup) {
+  int l1;
+  l1 = strchr(dup->structure, '&')-dup->structure;
+  printf("%s %3d,%-3d : %3d,%-3d (%5.2f)\n", dup->structure, dup->i+1-l1,
+         dup->i, dup->j, dup->j+strlen(dup->structure)-l1-2, dup->energy);
+}