RNAlib-2.1.1
data_structures.h
Go to the documentation of this file.
1 #ifndef __VIENNA_RNA_PACKAGE_DATA_STRUCTURES_H__
2 #define __VIENNA_RNA_PACKAGE_DATA_STRUCTURES_H__
3 
4 #include "energy_const.h"
10 /* to use floats instead of doubles in pf_fold() comment next line */
11 #define LARGE_PF
12 #ifdef LARGE_PF
13 #define FLT_OR_DBL double
14 #else
15 #define FLT_OR_DBL float
16 #endif
17 
18 #ifndef NBASES
19 #define NBASES 8
20 #endif
21 
22 #ifndef MAXALPHA
23 
26 #define MAXALPHA 20
27 #endif
28 
32 #define MAXDOS 1000
33 
34 #define VRNA_GQUAD_MAX_STACK_SIZE 7
35 #define VRNA_GQUAD_MIN_STACK_SIZE 2
36 #define VRNA_GQUAD_MAX_LINKER_LENGTH 15
37 #define VRNA_GQUAD_MIN_LINKER_LENGTH 1
38 #define VRNA_GQUAD_MIN_BOX_SIZE ((4*VRNA_GQUAD_MIN_STACK_SIZE)+(3*VRNA_GQUAD_MIN_LINKER_LENGTH))
39 #define VRNA_GQUAD_MAX_BOX_SIZE ((4*VRNA_GQUAD_MAX_STACK_SIZE)+(3*VRNA_GQUAD_MAX_LINKER_LENGTH))
40 
41 
42 /*
43 * ############################################################
44 * Here are the type definitions of various datastructures
45 * shared among the Vienna RNA Package
46 * ############################################################
47 */
48 
52 typedef struct plist {
53  int i;
54  int j;
55  float p;
56  int type;
57 } plist;
58 
62 typedef struct cpair {
63  int i,j,mfe;
64  float p, hue, sat;
65 } cpair;
66 
71 typedef struct {
72  float X; /* X coords */
73  float Y; /* Y coords */
74 } COORDINATE;
75 
79 typedef struct sect {
80  int i;
81  int j;
82  int ml;
83 } sect;
84 
88 typedef struct bondT {
89  unsigned int i;
90  unsigned int j;
91 } bondT;
92 
96 typedef struct bondTEn {
97  int i;
98  int j;
99  int energy;
100 } bondTEn;
101 
106 typedef struct{
107  int dangles;
114  int noLP;
115  int noGU;
117  int logML;
118  int circ;
119  int gquad;
121 
125 typedef struct{
126  int id;
127  int stack[NBPAIRS+1][NBPAIRS+1];
128  int hairpin[31];
129  int bulge[MAXLOOP+1];
130  int internal_loop[MAXLOOP+1];
131  int mismatchExt[NBPAIRS+1][5][5];
132  int mismatchI[NBPAIRS+1][5][5];
133  int mismatch1nI[NBPAIRS+1][5][5];
134  int mismatch23I[NBPAIRS+1][5][5];
135  int mismatchH[NBPAIRS+1][5][5];
136  int mismatchM[NBPAIRS+1][5][5];
137  int dangle5[NBPAIRS+1][5];
138  int dangle3[NBPAIRS+1][5];
139  int int11[NBPAIRS+1][NBPAIRS+1][5][5];
140  int int21[NBPAIRS+1][NBPAIRS+1][5][5][5];
141  int int22[NBPAIRS+1][NBPAIRS+1][5][5][5][5];
142  int ninio[5];
143  double lxc;
144  int MLbase;
145  int MLintern[NBPAIRS+1];
146  int MLclosing;
147  int TerminalAU;
148  int DuplexInit;
149  int Tetraloop_E[200];
150  char Tetraloops[1401];
151  int Triloop_E[40];
152  char Triloops[241];
153  int Hexaloop_E[40];
154  char Hexaloops[1801];
155  int TripleC;
156  int MultipleCA;
157  int MultipleCB;
158  int gquad [VRNA_GQUAD_MAX_STACK_SIZE + 1]
159  [3*VRNA_GQUAD_MAX_LINKER_LENGTH + 1];
160 
161  double temperature;
165 } paramT;
166 
170 typedef struct{
171  int id;
172  double expstack[NBPAIRS+1][NBPAIRS+1];
173  double exphairpin[31];
174  double expbulge[MAXLOOP+1];
175  double expinternal[MAXLOOP+1];
176  double expmismatchExt[NBPAIRS+1][5][5];
177  double expmismatchI[NBPAIRS+1][5][5];
178  double expmismatch23I[NBPAIRS+1][5][5];
179  double expmismatch1nI[NBPAIRS+1][5][5];
180  double expmismatchH[NBPAIRS+1][5][5];
181  double expmismatchM[NBPAIRS+1][5][5];
182  double expdangle5[NBPAIRS+1][5];
183  double expdangle3[NBPAIRS+1][5];
184  double expint11[NBPAIRS+1][NBPAIRS+1][5][5];
185  double expint21[NBPAIRS+1][NBPAIRS+1][5][5][5];
186  double expint22[NBPAIRS+1][NBPAIRS+1][5][5][5][5];
187  double expninio[5][MAXLOOP+1];
188  double lxc;
189  double expMLbase;
190  double expMLintern[NBPAIRS+1];
191  double expMLclosing;
192  double expTermAU;
193  double expDuplexInit;
194  double exptetra[40];
195  double exptri[40];
196  double exphex[40];
197  char Tetraloops[1401];
198  double expTriloop[40];
199  char Triloops[241];
200  char Hexaloops[1801];
201  double expTripleC;
202  double expMultipleCA;
203  double expMultipleCB;
204  double expgquad[VRNA_GQUAD_MAX_STACK_SIZE + 1]
205  [3*VRNA_GQUAD_MAX_LINKER_LENGTH + 1];
206 
207  double kT;
208  double pf_scale;
210  double temperature;
211  double alpha;
220 } pf_paramT;
221 
222 
223 
224 /*
225 * ############################################################
226 * SUBOPT data structures
227 * ############################################################
228 */
229 
230 
234 typedef struct {
235  int i;
236  int j;
237 } PAIR;
238 
242 typedef struct {
243  int i;
244  int j;
245  int array_flag;
246 } INTERVAL;
247 
251 typedef struct {
252  float energy;
253  char *structure;
254 } SOLUTION;
255 
256 /*
257 * ############################################################
258 * COFOLD data structures
259 * ############################################################
260 */
261 
265 typedef struct cofoldF {
266  /* free energies for: */
267  double F0AB;
268  double FAB;
269  double FcAB;
270  double FA;
271  double FB;
272 } cofoldF;
273 
277 typedef struct ConcEnt {
278  double A0;
279  double B0;
280  double ABc;
281  double AAc;
282  double BBc;
283  double Ac;
284  double Bc;
285 } ConcEnt;
286 
290 typedef struct pairpro{
291  struct plist *AB;
292  struct plist *AA;
293  struct plist *A;
294  struct plist *B;
295  struct plist *BB;
296 }pairpro;
297 
308 typedef struct {
309  unsigned i;
310  unsigned j;
311  float p;
312  float ent;
313  short bp[8];
314  char comp;
315 } pair_info;
316 
317 
318 /*
319 * ############################################################
320 * FINDPATH data structures
321 * ############################################################
322 */
323 
327 typedef struct move {
328  int i; /* i,j>0 insert; i,j<0 delete */
329  int j;
330  int when; /* 0 if still available, else resulting distance from start */
331  int E;
332 } move_t;
333 
337 typedef struct intermediate {
338  short *pt;
339  int Sen;
340  int curr_en;
343 
347 typedef struct path {
348  double en;
349  char *s;
350 } path_t;
351 
352 /*
353 * ############################################################
354 * RNAup data structures
355 * ############################################################
356 */
357 
361 typedef struct pu_contrib {
362  double **H;
363  double **I;
364  double **M;
365  double **E;
366  int length;
367  int w;
368 } pu_contrib;
369 
373 typedef struct interact {
374  double *Pi;
375  double *Gi;
376  double Gikjl;
378  double Gikjl_wo;
379  int i;
380  int k;
381  int j;
382  int l;
383  int length;
384 } interact;
385 
389 typedef struct pu_out {
390  int len;
391  int u_vals;
392  int contribs;
393  char **header;
394  double **u_values;
395 } pu_out;
396 
400 typedef struct constrain{
401  int *indx;
402  char *ptype;
403 } constrain;
404 
405 /*
406 * ############################################################
407 * RNAduplex data structures
408 * ############################################################
409 */
410 
414 typedef struct {
415  int i;
416  int j;
417  int end;
418  char *structure;
419  double energy;
420  double energy_backtrack;
421  double opening_backtrack_x;
422  double opening_backtrack_y;
423  int offset;
424  double dG1;
425  double dG2;
426  double ddG;
427  int tb;
428  int te;
429  int qb;
430  int qe;
431 } duplexT;
432 
433 /*
434 * ############################################################
435 * RNAsnoop data structures
436 * ############################################################
437 */
438 
442 typedef struct node {
443  int k;
444  int energy;
445  struct node *next;
446 } folden;
447 
451 typedef struct {
452  int i;
453  int j;
454  int u;
455  char *structure;
456  float energy;
457  float Duplex_El;
458  float Duplex_Er;
459  float Loop_E;
460  float Loop_D;
461  float pscd;
462  float psct;
463  float pscg;
464  float Duplex_Ol;
465  float Duplex_Or;
466  float Duplex_Ot;
467  float fullStemEnergy;
468 } snoopT;
469 
470 
471 
472 
473 
474 
475 
476 /*
477 * ############################################################
478 * PKplex data structures
479 * ############################################################
480 */
481 
485 typedef struct dupVar{
486  int i;
487  int j;
488  int end;
489  char *pk_helix;
490  char *structure;
491  double energy;
492  int offset;
493  double dG1;
494  double dG2;
495  double ddG;
496  int tb;
497  int te;
498  int qb;
499  int qe;
500  int inactive;
501  int processed;
502 } dupVar;
503 
504 
505 
506 /*
507 * ############################################################
508 * 2Dfold data structures
509 * ############################################################
510 */
511 
526 typedef struct{
527  int k;
528  int l;
529  float en;
530  char *s;
532 
538 typedef struct{
541  char *ptype;
542  char *sequence;
543  short *S, *S1;
544  unsigned int maxD1;
545  unsigned int maxD2;
548  unsigned int *mm1;
549  unsigned int *mm2;
551  int *my_iindx;
553  double temperature;
554 
555  unsigned int *referenceBPs1;
556  unsigned int *referenceBPs2;
557  unsigned int *bpdist;
559  short *reference_pt1;
560  short *reference_pt2;
561  int circ;
562  int dangles;
563  unsigned int seq_length;
564 
565  int ***E_F5;
566  int ***E_F3;
567  int ***E_C;
568  int ***E_M;
569  int ***E_M1;
570  int ***E_M2;
571 
572  int **E_Fc;
573  int **E_FcH;
574  int **E_FcI;
575  int **E_FcM;
576 
577  int **l_min_values;
578  int **l_max_values;
579  int *k_min_values;
580  int *k_max_values;
581 
582  int **l_min_values_m;
583  int **l_max_values_m;
584  int *k_min_values_m;
585  int *k_max_values_m;
586 
587  int **l_min_values_m1;
588  int **l_max_values_m1;
589  int *k_min_values_m1;
590  int *k_max_values_m1;
591 
592  int **l_min_values_f;
593  int **l_max_values_f;
594  int *k_min_values_f;
595  int *k_max_values_f;
596 
597  int **l_min_values_f3;
598  int **l_max_values_f3;
599  int *k_min_values_f3;
600  int *k_max_values_f3;
601 
602  int **l_min_values_m2;
603  int **l_max_values_m2;
604  int *k_min_values_m2;
605  int *k_max_values_m2;
606 
607  int *l_min_values_fc;
608  int *l_max_values_fc;
609  int k_min_values_fc;
610  int k_max_values_fc;
611 
612  int *l_min_values_fcH;
613  int *l_max_values_fcH;
614  int k_min_values_fcH;
615  int k_max_values_fcH;
616 
617  int *l_min_values_fcI;
618  int *l_max_values_fcI;
619  int k_min_values_fcI;
620  int k_max_values_fcI;
621 
622  int *l_min_values_fcM;
623  int *l_max_values_fcM;
624  int k_min_values_fcM;
625  int k_max_values_fcM;
626 
627  /* auxilary arrays for remaining set of coarse graining (k,l) > (k_max, l_max) */
628  int *E_F5_rem;
629  int *E_F3_rem;
630  int *E_C_rem;
631  int *E_M_rem;
632  int *E_M1_rem;
633  int *E_M2_rem;
634 
635  int E_Fc_rem;
636  int E_FcH_rem;
637  int E_FcI_rem;
638  int E_FcM_rem;
639 
640 #ifdef COUNT_STATES
641  unsigned long ***N_F5;
642  unsigned long ***N_C;
643  unsigned long ***N_M;
644  unsigned long ***N_M1;
645 #endif
646 } TwoDfold_vars;
647 
660 typedef struct{
661  int k;
662  int l;
663  FLT_OR_DBL q;
665 
672 typedef struct{
673 
674  unsigned int alloc;
675  char *ptype;
676  char *sequence;
677  short *S, *S1;
678  unsigned int maxD1;
679  unsigned int maxD2;
681  double temperature; /* temperature in last call to scale_pf_params */
682  double init_temp; /* temperature in last call to scale_pf_params */
683  FLT_OR_DBL *scale;
684  FLT_OR_DBL pf_scale;
685  pf_paramT *pf_params; /* holds all [unscaled] pf parameters */
686 
687  int *my_iindx;
688  int *jindx;
690  short *reference_pt1;
691  short *reference_pt2;
692 
693  unsigned int *referenceBPs1;
694  unsigned int *referenceBPs2;
695  unsigned int *bpdist;
697  unsigned int *mm1;
698  unsigned int *mm2;
700  int circ;
701  int dangles;
702  unsigned int seq_length;
703 
704  FLT_OR_DBL ***Q;
705  FLT_OR_DBL ***Q_B;
706  FLT_OR_DBL ***Q_M;
707  FLT_OR_DBL ***Q_M1;
708  FLT_OR_DBL ***Q_M2;
709 
710  FLT_OR_DBL **Q_c;
711  FLT_OR_DBL **Q_cH;
712  FLT_OR_DBL **Q_cI;
713  FLT_OR_DBL **Q_cM;
714 
715  int **l_min_values;
716  int **l_max_values;
717  int *k_min_values;
718  int *k_max_values;
719 
720  int **l_min_values_b;
721  int **l_max_values_b;
722  int *k_min_values_b;
723  int *k_max_values_b;
724 
725  int **l_min_values_m;
726  int **l_max_values_m;
727  int *k_min_values_m;
728  int *k_max_values_m;
729 
730  int **l_min_values_m1;
731  int **l_max_values_m1;
732  int *k_min_values_m1;
733  int *k_max_values_m1;
734 
735  int **l_min_values_m2;
736  int **l_max_values_m2;
737  int *k_min_values_m2;
738  int *k_max_values_m2;
739 
740  int *l_min_values_qc;
741  int *l_max_values_qc;
742  int k_min_values_qc;
743  int k_max_values_qc;
744 
745  int *l_min_values_qcH;
746  int *l_max_values_qcH;
747  int k_min_values_qcH;
748  int k_max_values_qcH;
749 
750  int *l_min_values_qcI;
751  int *l_max_values_qcI;
752  int k_min_values_qcI;
753  int k_max_values_qcI;
754 
755  int *l_min_values_qcM;
756  int *l_max_values_qcM;
757  int k_min_values_qcM;
758  int k_max_values_qcM;
759 
760  /* auxilary arrays for remaining set of coarse graining (k,l) > (k_max, l_max) */
761  FLT_OR_DBL *Q_rem;
762  FLT_OR_DBL *Q_B_rem;
763  FLT_OR_DBL *Q_M_rem;
764  FLT_OR_DBL *Q_M1_rem;
765  FLT_OR_DBL *Q_M2_rem;
766 
767  FLT_OR_DBL Q_c_rem;
768  FLT_OR_DBL Q_cH_rem;
769  FLT_OR_DBL Q_cI_rem;
770  FLT_OR_DBL Q_cM_rem;
771 
773 
774 #endif