1 /*
2 * PostScript and other output formats for RNA secondary structure plots
3 *
4 * c Ivo Hofacker, Peter F Stadler, Ronny Lorenz
5 * Vienna RNA package
6 */
7
8 #ifdef HAVE_CONFIG_H
9 #include "config.h"
10 #endif
11
12 #include <stdio.h>
13 #include <stdlib.h>
14 #include <math.h>
15 #include <string.h>
16 #include <ctype.h>
17 #include "ViennaRNA/utils/basic.h"
18 #include "ViennaRNA/utils/alignments.h"
19 #include "ViennaRNA/alphabet.h"
20 #include "ViennaRNA/model.h"
21 #include "ViennaRNA/fold_vars.h"
22 #include "ViennaRNA/gquad.h"
23 #include "ViennaRNA/plotting/layouts.h"
24 #include "ViennaRNA/plotting/utils.h"
25 #include "ViennaRNA/plotting/structures.h"
26 #include "ViennaRNA/plotting/RNApuzzler/RNApuzzler.h"
27 #include "ViennaRNA/plotting/RNApuzzler/RNAturtle.h"
28 #include "ViennaRNA/plotting/RNApuzzler/includes/svgArcs.inc"
29
30 #include "ViennaRNA/static/templates_postscript.h"
31 #include "ViennaRNA/static/templates_svg.h"
32
33 #include "ViennaRNA/plotting/ps_helpers.inc"
34 #include "ViennaRNA/plotting/svg_helpers.inc"
35
36 /*
37 #################################
38 # PRIVATE MACROS #
39 #################################
40 */
41 #define SIZE 452.
42
43 /*
44 #################################
45 # GLOBAL VARIABLES #
46 #################################
47 */
48
49 /*
50 #################################
51 # PRIVATE VARIABLES #
52 #################################
53 */
54
55 /*
56 #################################
57 # PRIVATE FUNCTION DECLARATIONS #
58 #################################
59 */
60 PRIVATE int
61 rnaplot_EPS(const char *seq,
62 const char *structure,
63 const char *ssfile,
64 const char *pre,
65 const char *post,
66 vrna_md_t *md_p,
67 vrna_plot_layout_t *layout);
68
69
70 /*
71 #################################
72 # BEGIN OF FUNCTION DEFINITIONS #
73 #################################
74 */
75 PUBLIC int
vrna_file_PS_rnaplot(const char * string,const char * structure,const char * ssfile,vrna_md_t * md_p)76 vrna_file_PS_rnaplot(const char *string,
77 const char *structure,
78 const char *ssfile,
79 vrna_md_t *md_p)
80 {
81 return vrna_file_PS_rnaplot_a(string, structure, ssfile, NULL, NULL, md_p);
82 }
83
84
85 PUBLIC int
vrna_file_PS_rnaplot_a(const char * seq,const char * structure,const char * ssfile,const char * pre,const char * post,vrna_md_t * md_p)86 vrna_file_PS_rnaplot_a(const char *seq,
87 const char *structure,
88 const char *ssfile,
89 const char *pre,
90 const char *post,
91 vrna_md_t *md_p)
92 {
93 int ret;
94 vrna_plot_layout_t *layout;
95
96 layout = vrna_plot_layout(structure, rna_plot_type);
97 ret = vrna_file_PS_rnaplot_layout(seq,
98 structure,
99 ssfile,
100 pre,
101 post,
102 md_p,
103 layout);
104
105 vrna_plot_layout_free(layout);
106
107 return ret;
108 }
109
110
111 PUBLIC int
vrna_file_PS_rnaplot_layout(const char * seq,const char * structure,const char * ssfile,const char * pre,const char * post,vrna_md_t * md_p,vrna_plot_layout_t * layout)112 vrna_file_PS_rnaplot_layout(const char *seq,
113 const char *structure,
114 const char *ssfile,
115 const char *pre,
116 const char *post,
117 vrna_md_t *md_p,
118 vrna_plot_layout_t *layout)
119 {
120 if (!ssfile) {
121 vrna_message_warning("vrna_file_PS_rnaplot*(): "
122 "Filename missing!");
123 } else if (!seq) {
124 vrna_message_warning("vrna_file_PS_rnaplot*(): "
125 "Sequence missing");
126 } else if (!structure) {
127 vrna_message_warning("vrna_file_PS_rnaplot*(): "
128 "Structure missing");
129 } else if (!layout) {
130 vrna_message_warning("vrna_file_PS_rnaplot*(): "
131 "Layout missing");
132 } else if ((strlen(seq) != strlen(structure)) ||
133 (strlen(structure) != layout->length)) {
134 vrna_message_warning("vrna_file_PS_rnaplot*(): "
135 "Sequence, structure, and coordinate layout have different lengths! "
136 "(%u vs. %u vs. %u)",
137 strlen(seq),
138 strlen(structure),
139 layout->length);
140 } else {
141 return rnaplot_EPS(seq,
142 structure,
143 ssfile,
144 pre,
145 post,
146 md_p,
147 layout);
148 }
149
150 return 0;
151 }
152
153
154 /*
155 #####################################
156 # BEGIN OF STATIC HELPER FUNCTIONS #
157 #####################################
158 */
159 PRIVATE int
rnaplot_EPS(const char * seq,const char * structure,const char * ssfile,const char * pre,const char * post,vrna_md_t * md_p,vrna_plot_layout_t * layout)160 rnaplot_EPS(const char *seq,
161 const char *structure,
162 const char *ssfile,
163 const char *pre,
164 const char *post,
165 vrna_md_t *md_p,
166 vrna_plot_layout_t *layout)
167 {
168 float xmin, xmax, ymin, ymax;
169 int i, length;
170 int ee, gb, ge, Lg, l[3], bbox[4];
171 float *X, *Y;
172 FILE *xyplot;
173 short *pair_table;
174 char *c, *string;
175 vrna_md_t md;
176
177 if (!md_p) {
178 set_model_details(&md);
179 md_p = &md;
180 }
181
182 string = strdup(seq);
183 length = strlen(string);
184
185 xyplot = fopen(ssfile, "w");
186 if (xyplot == NULL) {
187 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
188 return 0;
189 }
190
191 pair_table = vrna_ptable(structure);
192
193 bbox[0] = bbox[1] = 0;
194 bbox[2] = bbox[3] = 700;
195
196 print_PS_header(xyplot,
197 "RNA Secondary Structure Plot",
198 bbox,
199 md_p,
200 "To switch off outline pairs of sequence comment or\n"
201 "delete the appropriate line near the end of the file",
202 "RNAplot",
203 PS_MACRO_LAYOUT_BASE | ((pre || post) ? PS_MACRO_LAYOUT_EXTRAS : 0));
204
205 fprintf(xyplot, "%% data start here\n");
206
207 /* cut_point */
208 if ((c = strchr(structure, '&'))) {
209 int cutpoint;
210 cutpoint = c - structure;
211 string[cutpoint] = ' '; /* replace & with space */
212 fprintf(xyplot, "/cutpoint %d def\n", cutpoint);
213 }
214
215 /* sequence */
216 print_PS_sequence(xyplot, string);
217
218 /* coordinates */
219 print_PS_coords(xyplot, layout->x, layout->y, length);
220
221 fprintf(xyplot, "/arcs [\n");
222 if (layout->arcs) {
223 for (i = 0; i < length; i++) {
224 if (layout->arcs[6 * i + 2] > 0) {
225 /* 6*i+2 is the radius parameter */
226 fprintf(xyplot, "[%3.8f %3.8f %3.8f %3.8f %3.8f %3.8f]\n",
227 layout->arcs[6 * i + 0],
228 layout->arcs[6 * i + 1],
229 layout->arcs[6 * i + 2],
230 layout->arcs[6 * i + 3],
231 layout->arcs[6 * i + 4],
232 layout->arcs[6 * i + 5]
233 );
234 } else {
235 fprintf(xyplot, "[]\n");
236 }
237 }
238 } else {
239 for (i = 0; i < length; i++)
240 fprintf(xyplot, "[]\n");
241 }
242
243 fprintf(xyplot, "] def\n");
244
245 /* correction coordinates for quadratic beziers in case we produce a circplot */
246 if (rna_plot_type == VRNA_PLOT_TYPE_CIRCULAR)
247 fprintf(xyplot, "/cpr %6.2f def\n", (float)3 * length);
248
249 /* base pairs */
250 fprintf(xyplot, "/pairs [\n");
251 for (i = 1; i <= length; i++)
252 if (pair_table[i] > i)
253 fprintf(xyplot, "[%d %d]\n", i, pair_table[i]);
254
255 /* add gquad pairs */
256 ge = 0;
257 while ((ee = parse_gquad(structure + ge, &Lg, l)) > 0) {
258 int k;
259 fprintf(xyplot, "%% gquad\n");
260 ge += ee;
261 gb = ge - Lg * 4 - l[0] - l[1] - l[2] + 1; /* add pseudo-base pair encloding gquad */
262 for (k = 0; k < Lg; k++) {
263 int ii, jj, il;
264 for (il = 0, ii = gb + k; il < 3; il++) {
265 jj = ii + l[il] + Lg;
266 fprintf(xyplot, "[%d %d]\n", ii, jj);
267 ii = jj;
268 }
269 jj = gb + k;
270 fprintf(xyplot, "[%d %d]\n", jj, ii);
271 }
272 }
273
274 fprintf(xyplot, "] def\n\n");
275
276 fprintf(xyplot, "init\n\n");
277 /* draw the data */
278 if (pre) {
279 fprintf(xyplot, "%% Start Annotations\n");
280 fprintf(xyplot, "%s\n", pre);
281 fprintf(xyplot, "%% End Annotations\n");
282 }
283
284 fprintf(xyplot,
285 "%% switch off outline pairs or bases by removing these lines\n"
286 "drawoutline\n"
287 "drawpairs\n"
288 "drawbases\n");
289
290 if (post) {
291 fprintf(xyplot, "%% Start Annotations\n");
292 fprintf(xyplot, "%s\n", post);
293 fprintf(xyplot, "%% End Annotations\n");
294 }
295
296 print_PS_footer(xyplot);
297
298 fclose(xyplot);
299
300 free(string);
301 free(pair_table);
302
303 return 1; /* success */
304 }
305
306
307 /* options for gml output:
308 * uppercase letters: print sequence labels
309 * lowercase letters: no sequence lables
310 * graphics information:
311 * x X simple xy plot
312 * (nothing else implemented at present)
313 * default: no graphics data at all
314 */
315 PUBLIC int
gmlRNA(char * string,char * structure,char * ssfile,char option)316 gmlRNA(char *string,
317 char *structure,
318 char *ssfile,
319 char option)
320 {
321 FILE *gmlfile;
322 int i;
323 int length;
324 short *pair_table;
325 float *X, *Y;
326
327 gmlfile = fopen(ssfile, "w");
328 if (gmlfile == NULL) {
329 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
330 return 0;
331 }
332
333 length = strlen(string);
334
335 pair_table = vrna_ptable(structure);
336
337 switch (option) {
338 case 'X':
339 case 'x':
340 /* Simple XY Plot */
341 if (rna_plot_type == 0)
342 i = vrna_plot_coords_simple_pt(pair_table, &X, &Y);
343 else
344 i = vrna_plot_coords_naview_pt(pair_table, &X, &Y);
345
346 if (i != length)
347 vrna_message_warning("strange things happening in gmlRNA ...");
348
349 break;
350 default:
351 /* No Graphics Information */
352 X = NULL;
353 Y = NULL;
354 }
355
356 fprintf(gmlfile,
357 "# Vienna RNA Package %s\n"
358 "# GML Output\n"
359 "# CreationDate: %s\n"
360 "# Name: %s\n"
361 "# Options: %s\n", VERSION, vrna_time_stamp(), ssfile, option_string());
362 fprintf(gmlfile,
363 "graph [\n"
364 " directed 0\n");
365 for (i = 1; i <= length; i++) {
366 fprintf(gmlfile,
367 " node [ id %d ", i);
368 if (option)
369 fprintf(gmlfile,
370 "label \"%c\"", string[i - 1]);
371
372 if ((option == 'X') || (option == 'x'))
373 fprintf(gmlfile,
374 "\n graphics [ x %9.4f y %9.4f ]\n", X[i - 1], Y[i - 1]);
375
376 fprintf(gmlfile, " ]\n");
377 }
378 for (i = 1; i < length; i++)
379 fprintf(gmlfile,
380 "edge [ source %d target %d ]\n", i, i + 1);
381 for (i = 1; i <= length; i++) {
382 if (pair_table[i] > i)
383 fprintf(gmlfile,
384 "edge [ source %d target %d ]\n", i, pair_table[i]);
385 }
386 fprintf(gmlfile, "]\n");
387 fclose(gmlfile);
388
389 free(pair_table);
390 free(X);
391 free(Y);
392 return 1; /* success */
393 }
394
395
396 PUBLIC int
PS_rna_plot_snoop_a(const char * string,const char * structure,const char * ssfile,int * relative_access,const char * seqs[])397 PS_rna_plot_snoop_a(const char *string,
398 const char *structure,
399 const char *ssfile,
400 int *relative_access,
401 const char *seqs[])
402 {
403 int i, length, bbox[4];
404 float *X, *Y;
405 FILE *xyplot;
406 short *pair_table;
407 short *pair_table_snoop;
408 vrna_md_t md;
409
410 set_model_details(&md);
411
412 length = strlen(string);
413
414 xyplot = fopen(ssfile, "w");
415 if (xyplot == NULL) {
416 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
417 return 0;
418 }
419
420 pair_table = vrna_ptable(structure);
421 pair_table_snoop = vrna_pt_snoop_get(structure);
422
423 if (rna_plot_type == 0)
424 i = vrna_plot_coords_simple_pt(pair_table, &X, &Y);
425 else
426 i = vrna_plot_coords_naview_pt(pair_table, &X, &Y);
427
428 if (i != length)
429 vrna_message_warning("strange things happening in PS_rna_plot...");
430
431 /* printf("cut_point %d\n", cut_point); */
432
433 /*
434 * for (i = 1; i < length; i++) {
435 * printf("%d X %f Y %f \n", i, X[i], Y[i]);
436 * xmin = X[i] < xmin ? X[i] : xmin;
437 * xmax = X[i] > xmax ? X[i] : xmax;
438 * ymin = Y[i] < ymin ? Y[i] : ymin;
439 * ymax = Y[i] > ymax ? Y[i] : ymax;
440 * }
441 * localize centre of the interaction bucket. Geometry
442 */
443
444 for (i = 1; i < cut_point; i++) {
445 /* interior loop of size 0 */
446 if (pair_table_snoop[i] != 0) {
447 X[i - 1] = X[pair_table_snoop[i] - 1];
448 Y[i - 1] = Y[pair_table_snoop[i] - 1];
449 } else if (pair_table_snoop[i - 1] && pair_table_snoop[i + 1]) {
450 /* interior loop of size 1 */
451 X[i - 1] = X[pair_table_snoop[i - 1] - 1 - 1];
452 Y[i - 1] = Y[pair_table_snoop[i - 1] - 1 - 1];
453 } else if (pair_table_snoop[i - 1] && pair_table_snoop[i + 2]) {
454 /* interior loop of size 2 */
455 if (pair_table_snoop[i - 1] - pair_table_snoop[i + 2] == 2) {
456 X[i - 1] = X[pair_table_snoop[i - 1] - 2];
457 Y[i - 1] = Y[pair_table_snoop[i - 1] - 2];
458 X[i] = X[pair_table_snoop[i + 2]];
459 Y[i] = Y[pair_table_snoop[i + 2]];
460 i++;
461 } else if (pair_table[pair_table_snoop[i - 1] - 1]) {
462 X[i - 1] = X[pair_table_snoop[i - 1] - 2];
463 Y[i - 1] = Y[pair_table_snoop[i - 1] - 2];
464 X[i] = X[pair_table[pair_table_snoop[i - 1] - 1] - 1];
465 Y[i] = Y[pair_table[pair_table_snoop[i - 1] - 1] - 1];
466 i++;
467 } else if (pair_table[pair_table_snoop[i - 1] - 2]) {
468 X[i - 1] = X[pair_table_snoop[i - 1] - 3];
469 Y[i - 1] = Y[pair_table_snoop[i - 1] - 3];
470 X[i] = X[pair_table[pair_table_snoop[i - 1] - 2] - 1];
471 Y[i] = Y[pair_table[pair_table_snoop[i - 1] - 2] - 1];
472 i++;
473 } else if (pair_table[pair_table_snoop[i - 1] - 3]) {
474 X[i - 1] = X[pair_table_snoop[i - 1] - 4];
475 Y[i - 1] = Y[pair_table_snoop[i - 1] - 4];
476 X[i] = X[pair_table[pair_table_snoop[i - 1] - 3] - 1];
477 Y[i] = Y[pair_table[pair_table_snoop[i - 1] - 3] - 1];
478 i++;
479 } else {
480 X[i - 1] = X[pair_table_snoop[i - 1] - 2];
481 Y[i - 1] = Y[pair_table_snoop[i - 1] - 2];
482 X[i] = X[pair_table_snoop[i + 2]];
483 Y[i] = Y[pair_table_snoop[i + 2]];
484 i++;
485 }
486 } else if (pair_table_snoop[i - 1] && pair_table_snoop[i + 3]) {
487 /* interior loop of size 2 */
488 if (pair_table[pair_table_snoop[i - 1] - 1]) {
489 X[i - 1] = 0.5 * (X[pair_table_snoop[i - 1] - 1] + X[pair_table_snoop[i - 1] - 2]);
490 Y[i - 1] = 0.5 * (Y[pair_table_snoop[i - 1] - 1] + Y[pair_table_snoop[i - 1] - 2]);
491 X[i] = 0.5 *
492 (X[pair_table[pair_table_snoop[i - 1] - 1] - 1] +
493 X[pair_table_snoop[i - 1] - 2]);
494 Y[i] = 0.5 *
495 (Y[pair_table[pair_table_snoop[i - 1] - 1] - 1] + Y[pair_table_snoop[i - 1] - 2]);
496 X[i + 1] = 0.5 *
497 (X[pair_table[pair_table_snoop[i - 1] - 1] - 2] +
498 X[pair_table[pair_table_snoop[i - 1] - 1] - 1]);
499 Y[i + 1] = 0.5 *
500 (Y[pair_table[pair_table_snoop[i - 1] - 1] - 2] +
501 Y[pair_table[pair_table_snoop[i - 1] - 1] - 1]);
502 i++;
503 i++;
504 } else if (pair_table[pair_table_snoop[i - 1] - 2]) {
505 X[i - 1] = 0.5 * (X[pair_table_snoop[i - 1] - 2] + X[pair_table_snoop[i - 1] - 3]);
506 Y[i - 1] = 0.5 * (Y[pair_table_snoop[i - 1] - 2] + Y[pair_table_snoop[i - 1] - 3]);
507 X[i] = 0.5 *
508 (X[pair_table[pair_table_snoop[i - 1] - 2] - 1] +
509 X[pair_table_snoop[i - 1] - 3]);
510 Y[i] = 0.5 *
511 (Y[pair_table[pair_table_snoop[i - 1] - 2] - 1] + Y[pair_table_snoop[i - 1] - 3]);
512 X[i + 1] = 0.5 *
513 (X[pair_table[pair_table_snoop[i - 1] - 2] - 2] +
514 X[pair_table[pair_table_snoop[i - 1] - 2] - 1]);
515 Y[i + 1] = 0.5 *
516 (Y[pair_table[pair_table_snoop[i - 1] - 2] - 2] +
517 Y[pair_table[pair_table_snoop[i - 1] - 2] - 1]);
518 i++;
519 i++;
520 } else if (pair_table[pair_table_snoop[i - 1] - 3]) {
521 X[i - 1] = 0.5 * (X[pair_table_snoop[i - 1] - 3] + X[pair_table_snoop[i - 1] - 4]);
522 Y[i - 1] = 0.5 * (Y[pair_table_snoop[i - 1] - 3] + Y[pair_table_snoop[i - 1] - 4]);
523 X[i] = 0.5 *
524 (X[pair_table[pair_table_snoop[i - 1] - 3] - 1] +
525 X[pair_table_snoop[i - 1] - 4]);
526 Y[i] = 0.5 *
527 (Y[pair_table[pair_table_snoop[i - 1] - 3] - 1] + Y[pair_table_snoop[i - 1] - 4]);
528 X[i + 1] = 0.5 *
529 (X[pair_table[pair_table_snoop[i - 1] - 3] - 2] +
530 X[pair_table[pair_table_snoop[i - 1] - 3] - 1]);
531 Y[i + 1] = 0.5 *
532 (Y[pair_table[pair_table_snoop[i - 1] - 3] - 2] +
533 Y[pair_table[pair_table_snoop[i - 1] - 3] - 1]);
534 i++;
535 i++;
536 } else {
537 X[i - 1] = X[pair_table_snoop[i - 1] - 2];
538 Y[i - 1] = Y[pair_table_snoop[i - 1] - 2];
539 X[i] = X[pair_table_snoop[i - 1] - 2];
540 Y[i] = Y[pair_table_snoop[i - 1] - 2];
541 X[i + 1] = X[pair_table_snoop[i - 1] - 2];
542 Y[i + 1] = Y[pair_table_snoop[i - 1] - 2];
543 i++;
544 i++;
545 }
546 }
547 }
548 double xC;
549 double yC;
550 float X0 = -1, Y0 = -1, X1 = -1, Y1 = -1, X2 = -1, Y2 = -1;
551
552 /* int c1,c2,c3; */
553 for (i = 1; i < cut_point; i++) {
554 if (pair_table_snoop[i]) {
555 X0 = X[pair_table_snoop[i] - 1];
556 Y0 = Y[pair_table_snoop[i] - 1];
557 /* c1=pair_table_snoop[i]; */
558 i++;
559 break;
560 }
561 }
562 for (; i < cut_point; i++) {
563 if (pair_table_snoop[i]) {
564 X1 = X[pair_table_snoop[i] - 1];
565 Y1 = Y[pair_table_snoop[i] - 1];
566 /* c2=pair_table_snoop[i]; */
567 i++;
568 break;
569 }
570 }
571 for (; i < cut_point; i++) {
572 if (pair_table_snoop[i]) {
573 X2 = X[pair_table_snoop[i] - 1];
574 Y2 = Y[pair_table_snoop[i] - 1];
575 /* c3=pair_table_snoop[i]; */
576 i++;
577 break;
578 }
579 }
580 /*
581 * for(i=cut_point-2;i>pair_table_snoop[c1]; i--){
582 * if(pair_table_snoop[i]){
583 * X1=X[pair_table_snoop[i]-1];Y1=Y[pair_table_snoop[i]-1];
584 * c2=pair_table_snoop[i];
585 * i++;
586 * break;
587 * }
588 * }
589 * for(i=pair_table_snoop[c1]+1;i<pair_table_snoop[c2]; i++){
590 * if(pair_table_snoop[i]){
591 * X2=X[pair_table_snoop[i]-1];Y2=Y[pair_table_snoop[i]-1];
592 * c3=pair_table_snoop[i];
593 * i++;
594 * break;
595 * }
596 * }
597 */
598 if (X0 < 0 || X1 < 0 || X2 < 0) {
599 printf("Could not get the center of the binding bucket. No ps file will be produced!\n");
600 fclose(xyplot);
601 free(pair_table);
602 free(pair_table_snoop);
603 free(X);
604 free(Y);
605 pair_table = NULL;
606 pair_table_snoop = NULL;
607 X = NULL;
608 Y = NULL;
609 return 0;
610 }
611
612 double alpha = (X0 - X1) / (Y1 - Y0);
613 double alpha_p = (X1 - X2) / (Y2 - Y1);
614 double b = (Y0 + Y1 - alpha * (X0 + X1)) * 0.5;
615 double b_p = (Y1 + Y2 - alpha_p * (X1 + X2)) * 0.5;
616
617 /* if(abs(alpha -alpha_p) > 0.0000001){ */
618 xC = (b_p - b) / (alpha - alpha_p);
619 yC = alpha * xC + b;
620 for (i = 1; i < cut_point; i++) {
621 X[i - 1] = X[i - 1] + 0.25 * (xC - X[i - 1]);
622 Y[i - 1] = Y[i - 1] + 0.25 * (yC - Y[i - 1]);
623 }
624
625 bbox[0] = bbox[1] = 0;
626 bbox[2] = bbox[3] = 700;
627
628 print_PS_header(xyplot,
629 "RNA Secondary Structure Plot",
630 bbox,
631 &md,
632 "To switch off outline pairs of sequence comment or\n"
633 "delete the appropriate line near the end of the file",
634 "RNAplot",
635 PS_MACRO_LAYOUT_BASE | PS_MACRO_LAYOUT_EXTRAS);
636
637 char **A;
638
639 if (seqs) {
640 A = vrna_annotate_covar_db_extended((const char **)seqs,
641 structure,
642 &md,
643 VRNA_BRACKETS_RND | VRNA_BRACKETS_ANG | VRNA_BRACKETS_SQR);
644 }
645
646 fprintf(xyplot, "%% data start here\n");
647 /* cut_point */
648 if (cut_point > 0 && cut_point <= strlen(string))
649 fprintf(xyplot, "/cutpoint %d def\n", cut_point - 1);
650
651 /* sequence */
652 print_PS_sequence(xyplot, string);
653
654 /* coordinates */
655 print_PS_coords(xyplot, X, Y, length);
656
657 /* base pairs */
658 fprintf(xyplot, "/pairs [\n");
659 for (i = 1; i <= length; i++)
660 if (pair_table[i] > i)
661 fprintf(xyplot, "[%d %d]\n", i, pair_table[i]);
662
663 for (i = 1; i <= length; i++)
664 if (pair_table_snoop[i] > i)
665 fprintf(xyplot, "[%d %d]\n", i, pair_table_snoop[i]);
666
667 fprintf(xyplot, "] def\n\n");
668 if (relative_access) {
669 fprintf(xyplot, "/S [\n");
670 for (i = 0; i < cut_point - 1; i++)
671 fprintf(xyplot, " %f\n", (float)relative_access[i] / 100);
672 fprintf(xyplot, "]\n bind def\n");
673 fprintf(xyplot, "/invert false def\n");
674 fprintf(xyplot, "/range 0.8 def\n");
675 fprintf(xyplot, "/drawreliability {\n"
676 "/Smax 2.6 def\n"
677 " 0 \n"
678 " coor 0 cutpoint getinterval {\n"
679 " aload pop\n"
680 " S 3 index get\n"
681 " Smax div range mul\n"
682 " invert {range exch sub} if\n"
683 " 1 1 sethsbcolor\n"
684 " newpath\n"
685 " fsize 2.5 div 0 360 arc\n"
686 " fill\n"
687 " 1 add\n"
688 " } forall\n"
689 "\n"
690 "} bind def\n");
691 }
692
693 fprintf(xyplot, "init\n\n");
694 /*raw the data */
695 if (seqs) {
696 fprintf(xyplot, "%% Start Annotations\n");
697 fprintf(xyplot, "%s\n", A[0]);
698 fprintf(xyplot, "%% End Annotations\n");
699 }
700
701 fprintf(xyplot, "%%switch off outline pairs or bases by removing these lines\n");
702 if (relative_access)
703 fprintf(xyplot, "drawreliability\n");
704
705 fprintf(xyplot,
706 "drawoutline\n"
707 "drawpairs\n"
708 "drawbases\n");
709 /*
710 * fprintf(xyplot, "%d cmark\n",c1);
711 * fprintf(xyplot, "%d cmark\n",c2);
712 * fprintf(xyplot, "%d cmark\n",c3);
713 */
714 if (seqs) {
715 fprintf(xyplot, "%% Start Annotations\n");
716 fprintf(xyplot, "%s\n", A[1]);
717 fprintf(xyplot, "%% End Annotations\n");
718 }
719
720 print_PS_footer(xyplot);
721
722 fclose(xyplot);
723 if (seqs) {
724 free(A[0]);
725 free(A[1]);
726 free(A);
727 }
728
729 free(pair_table);
730 free(pair_table_snoop);
731 free(X);
732 free(Y);
733 return 1; /* success */
734 }
735
736
737 int
svg_rna_plot(char * string,char * structure,char * ssfile)738 svg_rna_plot(char *string,
739 char *structure,
740 char *ssfile)
741 {
742 float xmin, xmax, ymin, ymax, size;
743 int i, length;
744 float *X = NULL, *Y = NULL, *R = NULL, *CX = NULL, *CY = NULL;
745 FILE *xyplot;
746 short *pair_table;
747 double *arccoords = NULL, *arccoordsSVG = NULL;
748
749 length = strlen(string);
750
751 xyplot = fopen(ssfile, "w");
752 if (xyplot == NULL) {
753 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
754 return 0;
755 }
756
757 pair_table = vrna_ptable(structure);
758
759 switch (rna_plot_type) {
760 case VRNA_PLOT_TYPE_SIMPLE:
761 i = vrna_plot_coords_simple_pt(pair_table, &X, &Y);
762 break;
763 case VRNA_PLOT_TYPE_CIRCULAR:
764 {
765 int radius = 3 * length;
766 int dr = 0;
767 R = (float *)vrna_alloc((length + 1) * sizeof(float));
768 CX = (float *)vrna_alloc((length + 1) * sizeof(float));
769 CY = (float *)vrna_alloc((length + 1) * sizeof(float));
770 i = vrna_plot_coords_circular_pt(pair_table, &X, &Y);
771 for (i = 0; i < length; i++) {
772 if (i + 1 < pair_table[i + 1]) {
773 dr =
774 (pair_table[i + 1] - i + 1 <=
775 (length / 2 + 1)) ? pair_table[i + 1] - i : i + length - pair_table[i + 1];
776 R[i] = 1. - (2. * dr / (float)length);
777 } else if (pair_table[i + 1]) {
778 R[i] = R[pair_table[i + 1] - 1];
779 } else {
780 R[i] = 1.0;
781 }
782
783 CX[i] = X[i] * radius * R[i] + radius;
784 CY[i] = Y[i] * radius * R[i] + radius;
785 X[i] *= radius;
786 X[i] += radius;
787 Y[i] *= radius;
788 Y[i] += radius;
789 }
790 }
791 break;
792 case VRNA_PLOT_TYPE_TURTLE:
793 {
794 i = vrna_plot_coords_puzzler_pt(pair_table, &X, &Y, &arccoords, NULL);
795 transformPSArcsToSVG(i, arccoords, &arccoordsSVG); /* The arccoords for PS transformed to SVG Arcs data. Method implemented in ViennaRNA/plotting/RNApuzzler/includes/svgArcs.inc */
796 }
797 break;
798 case VRNA_PLOT_TYPE_PUZZLER:
799 {
800 i = vrna_plot_coords_puzzler_pt(pair_table, &X, &Y, &arccoords, NULL);
801 transformPSArcsToSVG(i, arccoords, &arccoordsSVG); /* See Turtle from above. */
802 }
803 break;
804
805 default:
806 i = vrna_plot_coords_naview_pt(pair_table, &X, &Y);
807 break;
808 }
809
810 if (i != length)
811 vrna_message_warning("strange things happening in PS_rna_plot...");
812
813 xmin = xmax = X[0];
814 ymin = ymax = Y[0];
815 for (i = 1; i < length; i++) {
816 xmin = X[i] < xmin ? X[i] : xmin;
817 xmax = X[i] > xmax ? X[i] : xmax;
818 ymin = Y[i] < ymin ? Y[i] : ymin;
819 ymax = Y[i] > ymax ? Y[i] : ymax;
820 }
821 for (i = 0; i < length; i++)
822 Y[i] = ymin + ymax - Y[i]; /* mirror coordinates so they look as in PS */
823
824 if (rna_plot_type == VRNA_PLOT_TYPE_CIRCULAR)
825 for (i = 0; i < length; i++) {
826 CY[i] = ymin + ymax - CY[i]; /* mirror coordinates so they look as in PS */
827 }
828
829 size = MAX2((xmax - xmin), (ymax - ymin));
830 size += 15; /* add some so the bounding box isn't too tight */
831
832 float scale[2] = {
833 SIZE / size, SIZE / size
834 };
835 float transform[2] = {
836 (float)(size - xmin - xmax) / 2., (float)(size - ymin - ymax) / 2.
837 };
838
839 print_SVG_header(xyplot, scale, transform, 0);
840
841 if (!(rna_plot_type == VRNA_PLOT_TYPE_PUZZLER || rna_plot_type == VRNA_PLOT_TYPE_TURTLE)) {
842 print_SVG_backbone(xyplot, X, Y, length);
843 } else {
844 /*
845 * ################################### SVG for Puzzler and Turtle ##########################################
846 * Draw Backbone
847 * Idea: We look at a point and test if there exist an arc between previous point and the actual point.
848 * But beware: Most of the indexes of this part of the code are "half" magic numbers;they were made
849 * with brute force(trial and error) until it eventually worked.
850 */
851 short newLine = 0; /* Flag if new Polyline should be created */
852 fprintf(xyplot,
853 " <polyline class=\"backbone\" id=\"outline\" points=\"\n");
854 for (int j = 1; j <= length; j++) {
855 if (arccoordsSVG[2 * (j - 1)] < 0) {
856 /* No arc(no valid radius) -> draw backbone */
857 if (newLine) {
858 newLine = 0;
859 fprintf(xyplot,
860 " <polyline class=\"backbone\" id=\"outline%i\" points=\"\n",
861 j);
862 fprintf(xyplot, " %3.3f,%3.3f\n", X[j - 2], Y[j - 2]); /* If new line is created, it should include the previous previous point */
863 }
864
865 fprintf(xyplot, " %3.3f,%3.3f\n", X[j - 1], Y[j - 1]); /* Then we include the previous point */
866 } else {
867 /*
868 * If there exist an arc, and there is no newline, then we are the last point
869 * in the polyline -> we set newLine to True and end Polyline
870 */
871 if (!newLine) {
872 newLine = 1;
873 fprintf(xyplot, " \" />\n");
874 }
875 }
876 }
877 fprintf(xyplot, " \" />\n");
878
879
880 /* arcs */
881 fprintf(xyplot, " <g id=\"arcs\">\n");
882 for (int j = 0; j < length - 1; j++) {
883 if (arccoordsSVG[2 * (j + 1)] > 0) {
884 /* If arc exists, then we draw the arc. Not much more to say */
885 fprintf(xyplot,
886 " <path class=\"backbone\" d=\"M %6.5f, %6.5f A %6.5f,%6.5f, %6.5f,%i, %i, %6.5f, %6.5f\" />\n",
887 X[j],
888 Y[j],
889 arccoordsSVG[2 * (j + 1)],
890 arccoordsSVG[2 * (j + 1)],
891 0.,
892 0,
893 (int)arccoordsSVG[2 * (j + 1) + 1],
894 X[j + 1],
895 Y[j + 1]);
896 }
897 }
898 fprintf(xyplot, " </g>\n");
899
900 /* ################################# END SVG Puzzler and Turtle ######################################### */
901 }
902
903 print_SVG_pairs(xyplot, pair_table, X, Y, CX, CY, length, rna_plot_type);
904
905 print_SVG_bases(xyplot, X, Y, string, length);
906
907 print_SVG_footer(xyplot);
908
909 fclose(xyplot);
910
911 free(pair_table);
912 free(X);
913 free(Y);
914 free(R);
915 free(CX);
916 free(CY);
917 free(arccoords);
918 free(arccoordsSVG);
919
920 return 1; /* success */
921 }
922
923
924 /*--------------------------------------------------------------------------*/
925
926 PUBLIC int
ssv_rna_plot(char * string,char * structure,char * ssfile)927 ssv_rna_plot(char *string,
928 char *structure,
929 char *ssfile)
930 {
931 /* produce input for the SStructView java applet */
932 FILE *ssvfile;
933 int i, bp;
934 int length;
935 short *pair_table;
936 float *X, *Y;
937 float xmin, xmax, ymin, ymax;
938
939 ssvfile = fopen(ssfile, "w");
940 if (ssvfile == NULL) {
941 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
942 return 0;
943 }
944
945 length = strlen(string);
946 pair_table = vrna_ptable(structure);
947
948 /* make coordinates */
949 if (rna_plot_type == 0)
950 i = vrna_plot_coords_simple_pt(pair_table, &X, &Y);
951 else
952 i = vrna_plot_coords_naview_pt(pair_table, &X, &Y);
953
954 if (i != length)
955 vrna_message_warning("strange things happening in ssv_rna_plot...");
956
957 /* make coords nonegative */
958 xmin = xmax = X[0];
959 ymin = ymax = Y[0];
960 for (i = 1; i < length; i++) {
961 xmin = X[i] < xmin ? X[i] : xmin;
962 xmax = X[i] > xmax ? X[i] : xmax;
963 ymin = Y[i] < ymin ? Y[i] : ymin;
964 ymax = Y[i] > ymax ? Y[i] : ymax;
965 }
966 if (xmin < 1) {
967 for (i = 0; i <= length; i++)
968 X[i] -= xmin - 1;
969 xmin = 1;
970 }
971
972 if (ymin < 1) {
973 for (i = 0; i <= length; i++)
974 Y[i] -= ymin - 1;
975 ymin = 1;
976 }
977
978 #if 0
979 {
980 float size, xoff, yoff;
981 float JSIZE = 500; /* size of the java applet window */
982 /* rescale coordinates, center on square of size HSIZE */
983 size = MAX2((xmax - xmin), (ymax - ymin));
984 xoff = (size - xmax + xmin) / 2;
985 yoff = (size - ymax + ymin) / 2;
986 for (i = 0; i <= length; i++) {
987 X[i] = (X[i] - xmin + xoff) * (JSIZE - 10) / size + 5;
988 Y[i] = (Y[i] - ymin + yoff) * (JSIZE - 10) / size + 5;
989 }
990 }
991 #endif
992 /* */
993
994 fprintf(ssvfile,
995 "# Vienna RNA Package %s\n"
996 "# SStructView Output\n"
997 "# CreationDate: %s\n"
998 "# Name: %s\n"
999 "# Options: %s\n", VERSION, vrna_time_stamp(), ssfile, option_string());
1000 for (i = 1; i <= length; i++)
1001 fprintf(ssvfile, "BASE\t%d\t%c\t%d\t%d\n",
1002 i, string[i - 1], (int)(X[i - 1] + 0.5), (int)(Y[i - 1] + 0.5));
1003 for (bp = 1, i = 1; i <= length; i++)
1004 if (pair_table[i] > i)
1005 fprintf(ssvfile, "BASE-PAIR\tbp%d\t%d\t%d\n", bp++, i, pair_table[i]);
1006
1007 fclose(ssvfile);
1008
1009 free(pair_table);
1010 free(X);
1011 free(Y);
1012 return 1; /* success */
1013 }
1014
1015
1016 /*---------------------------------------------------------------------------*/
1017 PUBLIC int
xrna_plot(char * string,char * structure,char * ssfile)1018 xrna_plot(char *string,
1019 char *structure,
1020 char *ssfile)
1021 {
1022 /* produce input for XRNA RNA drawing program */
1023 FILE *ss_file;
1024 int i;
1025 int length;
1026 short *pair_table;
1027 float *X, *Y;
1028
1029 ss_file = fopen(ssfile, "w");
1030 if (ss_file == NULL) {
1031 vrna_message_warning("can't open file %s - not doing xy_plot", ssfile);
1032 return 0;
1033 }
1034
1035 length = strlen(string);
1036 pair_table = vrna_ptable(structure);
1037
1038 /* make coordinates */
1039 if (rna_plot_type == 0)
1040 i = vrna_plot_coords_simple_pt(pair_table, &X, &Y);
1041 else
1042 i = vrna_plot_coords_naview_pt(pair_table, &X, &Y);
1043
1044 if (i != length)
1045 vrna_message_warning("strange things happening in xrna_plot...");
1046
1047 fprintf(ss_file,
1048 "# Vienna RNA Package %s, XRNA output\n"
1049 "# CreationDate: %s\n"
1050 "# Options: %s\n", VERSION, vrna_time_stamp(), option_string());
1051 for (i = 1; i <= length; i++)
1052 /* XRNA likes to have coordinate mirrored, so we use (-X, Y) */
1053 fprintf(ss_file, "%d %c %6.2f %6.2f %d %d\n", i, string[i - 1],
1054 -X[i - 1], Y[i - 1], (pair_table[i]?1:0), pair_table[i]);
1055 fclose(ss_file);
1056
1057 free(pair_table);
1058 free(X);
1059 free(Y);
1060 return 1; /* success */
1061 }
1062
1063
1064 #ifndef VRNA_DISABLE_BACKWARD_COMPATIBILITY
1065
1066 /*
1067 * ###########################################
1068 * # deprecated functions below #
1069 *###########################################
1070 */
1071 PUBLIC int
PS_rna_plot(char * string,char * structure,char * ssfile)1072 PS_rna_plot(char *string,
1073 char *structure,
1074 char *ssfile)
1075 {
1076 return vrna_file_PS_rnaplot((const char *)string,
1077 (const char *)structure,
1078 (const char *)ssfile,
1079 NULL);
1080 }
1081
1082
1083 PUBLIC int
PS_rna_plot_a(char * string,char * structure,char * ssfile,char * pre,char * post)1084 PS_rna_plot_a(char *string,
1085 char *structure,
1086 char *ssfile,
1087 char *pre,
1088 char *post)
1089 {
1090 return vrna_file_PS_rnaplot_a((const char *)string,
1091 (const char *)structure,
1092 (const char *)ssfile,
1093 (const char *)pre,
1094 (const char *)post,
1095 NULL);
1096 }
1097
1098
1099 PUBLIC int
PS_rna_plot_a_gquad(char * string,char * structure,char * ssfile,char * pre,char * post)1100 PS_rna_plot_a_gquad(char *string,
1101 char *structure,
1102 char *ssfile,
1103 char *pre,
1104 char *post)
1105 {
1106 return vrna_file_PS_rnaplot_a((const char *)string,
1107 (const char *)structure,
1108 (const char *)ssfile,
1109 (const char *)pre,
1110 (const char *)post,
1111 NULL);
1112 }
1113
1114
1115 #endif
1116