binaries/src/clustalo/src/squid/revcomp.c

   1 /*****************************************************************
   2  * SQUID - a library of functions for biological sequence analysis
   3  * Copyright (C) 1992-2002 Washington University School of Medicine
   4  *
   5  *     This source code is freely distributed under the terms of the
   6  *     GNU General Public License. See the files COPYRIGHT and LICENSE
   7  *     for details.
   8  *****************************************************************/
   9
  10 /* revcomp.c
  11  *
  12  * Reverse complement of a IUPAC character string
  13  * RCS $Id: revcomp.c 217 2011-03-19 10:27:10Z andreas $ (Original squid RCS Id: revcomp.c,v 1.5 2002/06/25 20:06:06 eddy Exp)
  14  */
  15
  16 #include <stdio.h>
  17 #include <string.h>
  18 #include <ctype.h>
  19 #include "squid.h"
  20
  21 /* Function: revcomp()
  22  *
  23  * Purpose:  Reverse complement seq; store in comp.
  24  *           Can revcomp "in place" (revcomp(seq, seq)).
  25  *
  26  * Args:     comp  - destination for reverse complement of seq
  27  *           seq   - sequence to reverse complement
  28  *
  29  * Returns:  NULL on failure; or a (useless) pointer to comp.
  30  */
  31 char *
  32 revcomp(char *comp, char *seq)
  33 {
  34   char *s;
  35   char  c;
  36
  37   if (comp == NULL) return NULL;
  38   if (seq == NULL)  return NULL;
  39
  40   StrReverse(comp, seq);
  41   for (s = comp; *s != '\0'; s++)
  42     {
  43       c = *s;
  44       c = sre_toupper(c);
  45       switch (c) {
  46       case 'A': c = 'T'; break;
  47       case 'C': c = 'G'; break;
  48       case 'G': c = 'C'; break;
  49       case 'T': c = 'A'; break;
  50       case 'U': c = 'A'; break;
  51       case 'R': c = 'Y'; break;
  52       case 'Y': c = 'R'; break;
  53       case 'M': c = 'K'; break;
  54       case 'K': c = 'M'; break;
  55       case 'S': c = 'S'; break;
  56       case 'W': c = 'W'; break;
  57       case 'H': c = 'D'; break;
  58       case 'D': c = 'H'; break;
  59       case 'B': c = 'V'; break;
  60       case 'V': c = 'B'; break;
  61       default:  break;          /* anything else? leave it; it's prob a gap or an X */
  62       }
  63       if (islower((int) *s)) c = (char) sre_tolower((int) c);
  64       *s = c;
  65     }
  66   return comp;
  67 }
  68
  69 #ifdef REVCOMP_TESTDRIVER
  70 /* gcc -g -DREVCOMP_TESTDRIVER revcomp.c sre_string.c shuffle.c sre_math.c sre_ctype.c sqerror.c -lm
  71 */
  72 int
  73 main(void)
  74 {
  75   float p[4]     = {0.25, 0.25, 0.25, 0.25};
  76   char *alphabet = "ACGT";
  77   int   len      = 10;
  78   char *seq;
  79
  80   seq = RandomSequence(alphabet, p, 4, len);
  81   printf("%s\n", seq);
  82   revcomp(seq, seq);
  83   printf("%s\n", seq);
  84   free(seq);
  85   exit(0);
  86 }
  87 #endif