RNAlib-2.4.14
pair_mat.h
1 #ifndef VIENNA_RNA_PACKAGE_PAIR_MAT_H
2 #define VIENNA_RNA_PACKAGE_PAIR_MAT_H
3 
4 #include <ctype.h>
6 #include <ViennaRNA/fold_vars.h>
7 
8 #define NBASES 8
9 /*@notnull@*/
10 
11 #ifndef INLINE
12 # ifdef __GNUC__
13 # define INLINE inline
14 # else
15 # define INLINE
16 # endif
17 #endif
18 
19 static const char Law_and_Order[] = "_ACGUTXKI";
20 static int BP_pair[NBASES][NBASES] =
21  /* _ A C G U X K I */
22 { { 0, 0, 0, 0, 0, 0, 0, 0 },
23  { 0, 0, 0, 0, 5, 0, 0, 5 },
24  { 0, 0, 0, 1, 0, 0, 0, 0 },
25  { 0, 0, 2, 0, 3, 0, 0, 0 },
26  { 0, 6, 0, 4, 0, 0, 0, 6 },
27  { 0, 0, 0, 0, 0, 0, 2, 0 },
28  { 0, 0, 0, 0, 0, 1, 0, 0 },
29  { 0, 6, 0, 0, 5, 0, 0, 0 } };
30 
31 #define MAXALPHA 20 /* maximal length of alphabet */
32 
33 static short alias[MAXALPHA + 1];
34 static int pair[MAXALPHA + 1][MAXALPHA + 1];
35 /* rtype[pair[i][j]]:=pair[j][i] */
36 static int rtype[8] = {
37  0, 2, 1, 4, 3, 6, 5, 7
38 };
39 
40 #ifdef _OPENMP
41 #pragma omp threadprivate(Law_and_Order, BP_pair, alias, pair, rtype)
42 #endif
43 
44 /* for backward compatibility */
45 #define ENCODE(c) encode_char(c)
46 
47 static INLINE int
48 encode_char(char c)
49 {
50  /* return numerical representation of base used e.g. in pair[][] */
51  int code;
52 
53  c = toupper(c);
54 
55  if (energy_set > 0) {
56  code = (int)(c - 'A') + 1;
57  } else {
58  const char *pos;
59  pos = strchr(Law_and_Order, c);
60  if (pos == NULL)
61  code = 0;
62  else
63  code = (int)(pos - Law_and_Order);
64 
65  if (code > 5)
66  code = 0;
67 
68  if (code > 4)
69  code--; /* make T and U equivalent */
70  }
71 
72  return code;
73 }
74 
75 
76 /*@+boolint +charint@*/
77 /*@null@*/
78 extern char *nonstandards;
79 
80 static INLINE void
81 make_pair_matrix(void)
82 {
83  int i, j;
84 
85  if (energy_set == 0) {
86  for (i = 0; i < 5; i++)
87  alias[i] = (short)i;
88  alias[5] = 3; /* X <-> G */
89  alias[6] = 2; /* K <-> C */
90  alias[7] = 0; /* I <-> default base '@' */
91  for (i = 0; i < NBASES; i++)
92  for (j = 0; j < NBASES; j++)
93  pair[i][j] = BP_pair[i][j];
94  if (noGU)
95  pair[3][4] = pair[4][3] = 0;
96 
97  if (nonstandards != NULL) {
98  /* allow nonstandard bp's */
99  for (i = 0; i < (int)strlen(nonstandards); i += 2)
100  pair[encode_char(nonstandards[i])]
101  [encode_char(nonstandards[i + 1])] = 7;
102  }
103 
104  for (i = 0; i < NBASES; i++)
105  for (j = 0; j < NBASES; j++)
106  rtype[pair[i][j]] = pair[j][i];
107  } else {
108  for (i = 0; i <= MAXALPHA; i++)
109  for (j = 0; j <= MAXALPHA; j++)
110  pair[i][j] = 0;
111  if (energy_set == 1) {
112  for (i = 1; i < MAXALPHA; ) {
113  alias[i++] = 3; /* A <-> G */
114  alias[i++] = 2; /* B <-> C */
115  }
116  for (i = 1; i < MAXALPHA; i++) {
117  pair[i][i + 1] = 2; /* AB <-> GC */
118  i++;
119  pair[i][i - 1] = 1; /* BA <-> CG */
120  }
121  } else if (energy_set == 2) {
122  for (i = 1; i < MAXALPHA; ) {
123  alias[i++] = 1; /* A <-> A*/
124  alias[i++] = 4; /* B <-> U */
125  }
126  for (i = 1; i < MAXALPHA; i++) {
127  pair[i][i + 1] = 5; /* AB <-> AU */
128  i++;
129  pair[i][i - 1] = 6; /* BA <-> UA */
130  }
131  } else if (energy_set == 3) {
132  for (i = 1; i < MAXALPHA - 2; ) {
133  alias[i++] = 3; /* A <-> G */
134  alias[i++] = 2; /* B <-> C */
135  alias[i++] = 1; /* C <-> A */
136  alias[i++] = 4; /* D <-> U */
137  }
138  for (i = 1; i < MAXALPHA - 2; i++) {
139  pair[i][i + 1] = 2; /* AB <-> GC */
140  i++;
141  pair[i][i - 1] = 1; /* BA <-> CG */
142  i++;
143  pair[i][i + 1] = 5; /* CD <-> AU */
144  i++;
145  pair[i][i - 1] = 6; /* DC <-> UA */
146  }
147  } else {
148  vrna_message_error("What energy_set are YOU using??");
149  }
150 
151  for (i = 0; i <= MAXALPHA; i++)
152  for (j = 0; j <= MAXALPHA; j++)
153  rtype[pair[i][j]] = pair[j][i];
154  }
155 }
156 
157 
158 static INLINE short *
159 encode_sequence(const char *sequence,
160  short how)
161 {
162  unsigned int i, l = (unsigned int)strlen(sequence);
163  short *S = (short *)vrna_alloc(sizeof(short) * (l + 2));
164 
165  switch (how) {
166  /* standard encoding as always used for S */
167  case 0:
168  for (i = 1; i <= l; i++) /* make numerical encoding of sequence */
169  S[i] = (short)encode_char(sequence[i - 1]);
170  S[l + 1] = S[1];
171  S[0] = (short)l;
172  break;
173  /* encoding for mismatches of nostandard bases (normally used for S1) */
174  case 1:
175  for (i = 1; i <= l; i++)
176  S[i] = alias[(short)encode_char(sequence[i - 1])];
177  S[l + 1] = S[1];
178  S[0] = S[l];
179  break;
180  }
181 
182  return S;
183 }
184 
185 
186 #endif /* VIENNA_RNA_PACKAGE_PAIR_MAT_H */
void * vrna_alloc(unsigned size)
Allocate space safely.
int energy_set
0 = BP; 1=any with GC; 2=any with AU-parameter
int noGU
Global switch to forbid/allow GU base pairs at all.
char * nonstandards
contains allowed non standard base pairs
#define MAXALPHA
Maximal length of alphabet.
Definition: model.h:166
void vrna_message_error(const char *format,...)
Print an error message and die.
Here all all declarations of the global variables used throughout RNAlib.
General utility- and helper-functions used throughout the ViennaRNA Package.