2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
47 #include "gromacs/commandline/pargs.h"
49 #include "gromacs/fileio/futil.h"
50 #include "gmx_fatal.h"
53 #include "gromacs/fileio/matio.h"
56 #include "gromacs/fileio/trxio.h"
73 char tickfont[STRLEN];
88 char leglabel[STRLEN];
89 char leg2label[STRLEN];
99 /* MUST correspond to char *legend[] in main() */
101 elSel, elBoth, elFirst, elSecond, elNone, elNR
104 /* MUST correspond to char *combine[] in main() */
106 ecSel, ecHalves, ecAdd, ecSub, ecMult, ecDiv, ecNR
109 void get_params(const char *mpin, const char *mpout, t_psrec *psr)
111 static const char *gmx_bools[BOOL_NR+1] = { "no", "yes", NULL };
112 /* this must correspond to t_rgb *linecolors[] below */
113 static const char *colors[] = { "none", "black", "white", NULL };
119 wi = init_warning(FALSE, 0);
123 inp = read_inpfile(mpin, &ninp, wi);
129 ETYPE("black&white", psr->bw, gmx_bools);
130 RTYPE("linewidth", psr->linewidth, 1.0);
131 STYPE("titlefont", psr->titfont, "Helvetica");
132 RTYPE("titlefontsize", psr->titfontsize, 20.0);
133 ETYPE("legend", psr->legend, gmx_bools);
134 STYPE("legendfont", psr->legfont, psr->titfont);
135 STYPE("legendlabel", psr->leglabel, "");
136 STYPE("legend2label", psr->leg2label, psr->leglabel);
137 RTYPE("legendfontsize", psr->legfontsize, 14.0);
138 RTYPE("xbox", psr->xboxsize, 0.0);
139 RTYPE("ybox", psr->yboxsize, 0.0);
140 RTYPE("matrixspacing", psr->boxspacing, 20.0);
141 RTYPE("xoffset", psr->xoffs, 0.0);
142 RTYPE("yoffset", psr->yoffs, psr->xoffs);
143 RTYPE("boxlinewidth", psr->boxlinewidth, psr->linewidth);
144 RTYPE("ticklinewidth", psr->ticklinewidth, psr->linewidth);
145 RTYPE("zerolinewidth", psr->zerolinewidth, psr->ticklinewidth);
146 ETYPE("x-lineat0value", psr->X.lineatzero, colors);
147 RTYPE("x-major", psr->X.major, NOTSET);
148 RTYPE("x-minor", psr->X.minor, NOTSET);
149 RTYPE("x-firstmajor", psr->X.offset, 0.0);
150 ETYPE("x-majorat0", psr->X.first, gmx_bools);
151 RTYPE("x-majorticklen", psr->X.majorticklen, 8.0);
152 RTYPE("x-minorticklen", psr->X.minorticklen, 4.0);
153 STYPE("x-label", psr->X.label, "");
154 RTYPE("x-fontsize", psr->X.fontsize, 16.0);
155 STYPE("x-font", psr->X.font, psr->titfont);
156 RTYPE("x-tickfontsize", psr->X.tickfontsize, 10.0);
157 STYPE("x-tickfont", psr->X.tickfont, psr->X.font);
158 ETYPE("y-lineat0value", psr->Y.lineatzero, colors);
159 RTYPE("y-major", psr->Y.major, psr->X.major);
160 RTYPE("y-minor", psr->Y.minor, psr->X.minor);
161 RTYPE("y-firstmajor", psr->Y.offset, psr->X.offset);
162 ETYPE("y-majorat0", psr->Y.first, gmx_bools);
163 RTYPE("y-majorticklen", psr->Y.majorticklen, psr->X.majorticklen);
164 RTYPE("y-minorticklen", psr->Y.minorticklen, psr->X.minorticklen);
165 STYPE("y-label", psr->Y.label, psr->X.label);
166 RTYPE("y-fontsize", psr->Y.fontsize, psr->X.fontsize);
167 STYPE("y-font", psr->Y.font, psr->X.font);
168 RTYPE("y-tickfontsize", psr->Y.tickfontsize, psr->X.tickfontsize);
169 STYPE("y-tickfont", psr->Y.tickfont, psr->Y.font);
173 write_inpfile(mpout, ninp, inp, TRUE, wi);
176 done_warning(wi, FARGS);
179 t_rgb black = { 0, 0, 0 };
180 t_rgb white = { 1, 1, 1 };
181 t_rgb red = { 1, 0, 0 };
182 t_rgb blue = { 0, 0, 1 };
183 #define BLACK (&black)
184 /* this must correspond to *colors[] in get_params */
185 t_rgb *linecolors[] = { NULL, &black, &white, NULL };
187 gmx_bool diff_maps(int nmap1, t_mapping *map1, int nmap2, t_mapping *map2)
190 gmx_bool bDiff, bColDiff = FALSE;
199 for (i = 0; i < nmap1; i++)
201 if (!matelmt_cmp(map1[i].code, map2[i].code))
205 if (strcmp(map1[i].desc, map2[i].desc) != 0)
209 if ((map1[i].rgb.r != map2[i].rgb.r) ||
210 (map1[i].rgb.g != map2[i].rgb.g) ||
211 (map1[i].rgb.b != map2[i].rgb.b))
216 if (!bDiff && bColDiff)
218 fprintf(stderr, "Warning: two colormaps differ only in RGB value, using one colormap.\n");
225 void leg_discrete(t_psdata ps, real x0, real y0, char *label,
226 real fontsize, char *font, int nmap, t_mapping map[])
232 boxhh = fontsize+DDD;
235 ps_strfont(ps, font, fontsize);
236 yhh = y0+fontsize+3*DDD;
237 if ((int)strlen(label) > 0)
239 ps_ctext(ps, x0, yhh, label, eXLeft);
241 ps_moveto(ps, x0, y0);
242 for (i = 0; (i < nmap); i++)
245 ps_rgb(ps, &(map[i].rgb));
246 ps_fillbox(ps, DDD, DDD, DDD+fontsize, boxhh-DDD);
248 ps_box(ps, DDD, DDD, DDD+fontsize, boxhh-DDD);
249 ps_ctext(ps, boxhh+2*DDD, fontsize/3, map[i].desc, eXLeft);
251 ps_moverel(ps, DDD, -fontsize/3);
255 void leg_continuous(t_psdata ps, real x0, real x, real y0, char *label,
256 real fontsize, char *font,
257 int nmap, t_mapping map[],
262 real yhh, boxxh, boxyh;
269 boxxh = (real)x/(real)(nmap-mapoffset);
270 if (boxxh > fontsize)
276 xx0 = x0-((nmap-mapoffset)*boxxh)/2.0;
278 for (i = 0; (i < nmap-mapoffset); i++)
280 ps_rgb(ps, &(map[i+mapoffset].rgb));
281 ps_fillbox(ps, xx0+i*boxxh, y0, xx0+(i+1)*boxxh, y0+boxyh);
283 ps_strfont(ps, font, fontsize);
285 ps_box(ps, xx0, y0, xx0+(nmap-mapoffset)*boxxh, y0+boxyh);
287 yhh = y0+boxyh+3*DDD;
288 ps_ctext(ps, xx0+boxxh/2, yhh, map[0].desc, eXCenter);
289 if ((int)strlen(label) > 0)
291 ps_ctext(ps, x0, yhh, label, eXCenter);
293 ps_ctext(ps, xx0+((nmap-mapoffset)*boxxh)
294 - boxxh/2, yhh, map[nmap-1].desc, eXCenter);
297 void leg_bicontinuous(t_psdata ps, real x0, real x, real y0, char *label1,
298 char *label2, real fontsize, char *font,
299 int nmap1, t_mapping map1[], int nmap2, t_mapping map2[])
301 real xx1, xx2, x1, x2;
303 x1 = x/(nmap1+nmap2)*nmap1; /* width of legend 1 */
304 x2 = x/(nmap1+nmap2)*nmap2; /* width of legend 2 */
305 xx1 = x0-(x2/2.0)-fontsize; /* center of legend 1 */
306 xx2 = x0+(x1/2.0)+fontsize; /* center of legend 2 */
307 x1 -= fontsize/2; /* adjust width */
308 x2 -= fontsize/2; /* adjust width */
309 leg_continuous(ps, xx1, x1, y0, label1, fontsize, font, nmap1, map1, 0);
310 leg_continuous(ps, xx2, x2, y0, label2, fontsize, font, nmap2, map2, 0);
313 static real box_height(t_matrix *mat, t_psrec *psr)
315 return mat->ny*psr->yboxsize;
318 static real box_dh(t_psrec *psr)
320 return psr->boxspacing;
323 #define IS_ONCE (i == nmat-1)
324 static real box_dh_top(gmx_bool bOnce, t_psrec *psr)
328 if (psr->bTitle || (psr->bTitleOnce && bOnce) )
330 dh = 2*psr->titfontsize;
340 static gmx_bool box_do_all_x_maj_ticks(t_psrec *psr)
342 return (psr->boxspacing > (1.5*psr->X.majorticklen));
345 static gmx_bool box_do_all_x_min_ticks(t_psrec *psr)
347 return (psr->boxspacing > (1.5*psr->X.minorticklen));
350 static void draw_boxes(t_psdata ps, real x0, real y0, real w,
351 int nmat, t_matrix mat[], t_psrec *psr)
356 char **xtick, **ytick;
357 real xx, yy, dy, xx00, yy00, offset_x, offset_y;
358 int i, j, x, y, ntx, nty, strlength;
360 /* Only necessary when there will be no y-labels */
365 ps_linewidth(ps, psr->boxlinewidth);
367 for (i = 0; (i < nmat); i++)
369 dy = box_height(&(mat[i]), psr);
370 ps_box(ps, x0-1, yy00-1, x0+w+1, yy00+dy+1);
371 yy00 += dy+box_dh(psr)+box_dh_top(IS_ONCE, psr);
374 /* Draw the ticks on the axes */
375 ps_linewidth(ps, psr->ticklinewidth);
378 for (i = 0; (i < nmat); i++)
380 if (mat[i].flags & MAT_SPATIAL_X)
390 if (mat[i].flags & MAT_SPATIAL_Y)
401 for (j = 0; (j < ntx); j++)
403 sprintf(buf, "%g", mat[i].axis_x[j]);
404 xtick[j] = strdup(buf);
406 ps_strfont(ps, psr->X.tickfont, psr->X.tickfontsize);
407 for (x = 0; (x < ntx); x++)
409 xx = xx00 + (x + offset_x)*psr->xboxsize;
410 if ( ( bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.major) ||
411 (psr->X.first && (x == 0))) &&
412 ( (i == 0) || box_do_all_x_maj_ticks(psr) ) )
414 /* Longer tick marks */
415 ps_line (ps, xx, yy00, xx, yy00-psr->X.majorticklen);
416 /* Plot label on lowest graph only */
420 yy00-DDD-psr->X.majorticklen-psr->X.tickfontsize*0.8,
424 else if (bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.minor) &&
425 ( (i == 0) || box_do_all_x_min_ticks(psr) ) )
427 /* Shorter tick marks */
428 ps_line(ps, xx, yy00, xx, yy00-psr->X.minorticklen);
430 else if (bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.major) )
432 /* Even shorter marks, only each X.major */
433 ps_line(ps, xx, yy00, xx, yy00-(psr->boxspacing/2));
436 ps_strfont(ps, psr->Y.tickfont, psr->Y.tickfontsize);
438 for (j = 0; (j < nty); j++)
440 sprintf(buf, "%g", mat[i].axis_y[j]);
441 ytick[j] = strdup(buf);
444 for (y = 0; (y < nty); y++)
446 yy = yy00 + (y + offset_y)*psr->yboxsize;
447 if (bRmod(mat[i].axis_y[y], psr->Y.offset, psr->Y.major) ||
448 (psr->Y.first && (y == 0)))
451 strlength = max(strlength, (int)strlen(ytick[y]));
452 ps_line (ps, xx00, yy, xx00-psr->Y.majorticklen, yy);
453 ps_ctext(ps, xx00-psr->Y.majorticklen-DDD,
454 yy-psr->Y.tickfontsize/3.0, ytick[y], eXRight);
456 else if (bRmod(mat[i].axis_y[y], psr->Y.offset, psr->Y.minor) )
459 ps_line(ps, xx00, yy, xx00-psr->Y.minorticklen, yy);
465 /* Label on Y-axis */
466 if (!psr->bYonce || i == nmat/2)
468 if (strlen(psr->Y.label) > 0)
470 mylab = psr->Y.label;
474 mylab = mat[i].label_y;
476 if (strlen(mylab) > 0)
478 ps_strfont(ps, psr->Y.font, psr->Y.fontsize);
480 xxx = x0-psr->X.majorticklen-psr->X.tickfontsize*strlength-DDD;
481 ps_ctext(ps, yy00+box_height(&mat[i], psr)/2.0, 612.5-xxx,
487 yy00 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
489 /* Label on X-axis */
490 if (strlen(psr->X.label) > 0)
492 mylab = psr->X.label;
496 mylab = mat[0].label_x;
498 if (strlen(mylab) > 0)
500 ps_strfont(ps, psr->X.font, psr->X.fontsize);
501 ps_ctext(ps, x0+w/2, y0-DDD-psr->X.majorticklen-psr->X.tickfontsize*FUDGE-
502 psr->X.fontsize, mylab, eXCenter);
506 static void draw_zerolines(t_psdata out, real x0, real y0, real w,
507 int nmat, t_matrix mat[], t_psrec *psr)
509 real xx, yy, dy, xx00, yy00;
514 ps_linewidth(out, psr->zerolinewidth);
515 for (i = 0; (i < nmat); i++)
517 dy = box_height(&(mat[i]), psr);
518 /* mat[i].axis_x and _y were already set by draw_boxes */
519 if (psr->X.lineatzero)
521 ps_rgb(out, linecolors[psr->X.lineatzero]);
522 for (x = 0; (x < mat[i].nx); x++)
524 xx = xx00+(x+0.7)*psr->xboxsize;
525 /* draw lines whenever tick label almost zero (e.g. next trajectory) */
526 if (x != 0 && x < mat[i].nx-1 &&
527 abs(mat[i].axis_x[x]) <
528 0.1*abs(mat[i].axis_x[x+1]-mat[i].axis_x[x]) )
530 ps_line (out, xx, yy00, xx, yy00+dy+2);
534 if (psr->Y.lineatzero)
536 ps_rgb(out, linecolors[psr->Y.lineatzero]);
537 for (y = 0; (y < mat[i].ny); y++)
539 yy = yy00+(y+0.7)*psr->yboxsize;
540 /* draw lines whenever tick label almost zero (e.g. next trajectory) */
541 if (y != 0 && y < mat[i].ny-1 &&
542 abs(mat[i].axis_y[y]) <
543 0.1*abs(mat[i].axis_y[y+1]-mat[i].axis_y[y]) )
545 ps_line (out, xx00, yy, xx00+w+2, yy);
549 yy00 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
553 static void box_dim(int nmat, t_matrix mat[], t_matrix *mat2, t_psrec *psr,
554 int elegend, gmx_bool bFrame,
555 real *w, real *h, real *dw, real *dh)
558 real ww, hh, dww, dhh;
564 for (i = 0; (i < nmat); i++)
566 ww = max(ww, mat[i].nx*psr->xboxsize);
567 hh += box_height(&(mat[i]), psr);
568 maxytick = max(maxytick, mat[i].nx);
572 if (mat[0].label_y[0])
574 dww += 2.0*(psr->Y.fontsize+DDD);
576 if (psr->Y.major > 0)
578 dww += psr->Y.majorticklen + DDD +
579 psr->Y.tickfontsize*(log(maxytick)/log(10.0));
581 else if (psr->Y.minor > 0)
583 dww += psr->Y.minorticklen;
586 if (mat[0].label_x[0])
588 dhh += psr->X.fontsize+2*DDD;
590 if ( /* fool emacs auto-indent */
591 (elegend == elBoth && (mat[0].legend[0] || mat2[0].legend[0])) ||
592 (elegend == elFirst && mat[0].legend[0]) ||
593 (elegend == elSecond && mat2[0].legend[0]) )
595 dhh += 2*(psr->legfontsize*FUDGE+2*DDD);
599 dhh += psr->legfontsize*FUDGE+2*DDD;
601 if (psr->X.major > 0)
603 dhh += psr->X.tickfontsize*FUDGE+2*DDD+psr->X.majorticklen;
605 else if (psr->X.minor > 0)
607 dhh += psr->X.minorticklen;
610 hh += (nmat-1)*box_dh(psr);
611 hh += box_dh_top(TRUE, psr);
614 hh += (nmat-1)*box_dh_top(FALSE, psr);
623 int add_maps(t_mapping **newmap,
624 int nmap1, t_mapping map1[], int nmap2, t_mapping map2[])
626 static char mapper[] = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789!@#$%^&*()-_=+{}|;:',<.>/?";
631 nsymbols = strlen(mapper);
633 if (nmap > nsymbols*nsymbols)
635 gmx_fatal(FARGS, "Not enough symbols to merge the two colormaps\n");
637 printf("Combining colormaps of %d and %d elements into one of %d elements\n",
640 for (j = 0; j < nmap1; j++)
642 map[j].code.c1 = mapper[j % nsymbols];
645 map[j].code.c2 = mapper[j/nsymbols];
647 map[j].rgb.r = map1[j].rgb.r;
648 map[j].rgb.g = map1[j].rgb.g;
649 map[j].rgb.b = map1[j].rgb.b;
650 map[j].desc = map1[j].desc;
652 for (j = 0; j < nmap2; j++)
655 map[k].code.c1 = mapper[k % nsymbols];
658 map[k].code.c2 = mapper[k/nsymbols];
660 map[k].rgb.r = map2[j].rgb.r;
661 map[k].rgb.g = map2[j].rgb.g;
662 map[k].rgb.b = map2[j].rgb.b;
663 map[k].desc = map2[j].desc;
670 void xpm_mat(const char *outf, int nmat, t_matrix *mat, t_matrix *mat2,
671 gmx_bool bDiag, gmx_bool bFirstDiag)
675 int i, j, k, x, y, col;
677 t_mapping *map = NULL;
679 out = ffopen(outf, "w");
681 for (i = 0; i < nmat; i++)
683 if (!mat2 || !diff_maps(mat[i].nmap, mat[i].map, mat2[i].nmap, mat2[i].map))
685 write_xpm_m(out, mat[0]);
689 nmap = add_maps(&map, mat[i].nmap, mat[i].map, mat2[i].nmap, mat2[i].map);
690 for (x = 0; (x < mat[i].nx); x++)
692 for (y = 0; (y < mat[i].nx); y++)
694 if ((x < y) || ((x == y) && bFirstDiag)) /* upper left -> map1 */
696 col = mat[i].matrix[x][y];
698 else /* lower right -> map2 */
700 col = mat[i].nmap+mat[i].matrix[x][y];
702 if ((bDiag) || (x != y))
704 mat[i].matrix[x][y] = col;
708 mat[i].matrix[x][y] = 0;
715 if (mat2 && (strcmp(mat[i].title, mat2[i].title) != 0))
717 sprintf(mat[i].title+strlen(mat[i].title), " / %s", mat2[i].title);
719 if (mat2 && (strcmp(mat[i].legend, mat2[i].legend) != 0))
721 sprintf(mat[i].legend+strlen(mat[i].legend), " / %s", mat2[i].legend);
723 write_xpm_m(out, mat[i]);
729 static void tick_spacing(int n, real axis[], real offset, char axisnm,
730 real *major, real *minor)
733 gmx_bool bTryAgain, bFive;
734 int i, j, t, f = 0, ten;
736 real major_fact[NFACT] = {5, 4, 2, 1};
737 real minor_fact[NFACT] = {5, 4, 4, 5};
739 /* start with interval between 10 matrix points: */
740 space = max(10*axis[1]-axis[0], axis[min(10, n-1)]-axis[0]);
741 /* get power of 10 */
742 ten = (int)ceil(log(space)/log(10))-1;
744 for (t = ten+2; t > ten-3 && bTryAgain; t--)
746 for (f = 0; f < NFACT && bTryAgain; f++)
748 space = pow(10, t) * major_fact[f];
749 /* count how many ticks we would get: */
751 for (j = 0; j < n; j++)
753 if (bRmod(axis[j], offset, space) )
758 /* do we have a reasonable number of ticks ? */
759 bTryAgain = (i > min(10, n-1)) || (i < 5);
764 space = max(10*axis[1]-axis[0], axis[min(10, n-1)]-axis[0]);
765 fprintf(stderr, "Auto tick spacing failed for %c-axis, guessing %g\n",
769 *minor = space / minor_fact[f-1];
770 fprintf(stderr, "Auto tick spacing for %c-axis: major %g, minor %g\n",
771 axisnm, *major, *minor);
774 void ps_mat(const char *outf, int nmat, t_matrix mat[], t_matrix mat2[],
775 gmx_bool bFrame, gmx_bool bDiag, gmx_bool bFirstDiag,
776 gmx_bool bTitle, gmx_bool bTitleOnce, gmx_bool bYonce, int elegend,
777 real size, real boxx, real boxy, const char *m2p, const char *m2pout,
781 char buf[256], *legend;
785 int i, j, x, y, col, leg = 0;
788 int nmap1 = 0, nmap2 = 0, leg_nmap;
789 t_mapping *map1 = NULL, *map2 = NULL, *leg_map;
790 gmx_bool bMap1, bNextMap1, bDiscrete;
793 libm2p = m2p ? gmxlibfn(m2p) : m2p;
794 get_params(libm2p, m2pout, &psrec);
798 if (psr->X.major <= 0)
800 tick_spacing((mat[0].flags & MAT_SPATIAL_X) ? mat[0].nx + 1 : mat[0].nx,
801 mat[0].axis_x, psr->X.offset, 'X',
802 &(psr->X.major), &(psr->X.minor) );
804 if (psr->X.minor <= 0)
806 psr->X.minor = psr->X.major / 2;
808 if (psr->Y.major <= 0)
810 tick_spacing((mat[0].flags & MAT_SPATIAL_Y) ? mat[0].ny + 1 : mat[0].ny,
811 mat[0].axis_y, psr->Y.offset, 'Y',
812 &(psr->Y.major), &(psr->Y.minor) );
814 if (psr->Y.minor <= 0)
816 psr->Y.minor = psr->Y.major / 2;
821 psr->xboxsize = boxx;
822 psr->yboxsize = boxx;
826 psr->yboxsize = boxy;
829 if (psr->xboxsize == 0)
831 psr->xboxsize = size/mat[0].nx;
832 printf("Set the x-size of the box to %.3f\n", psr->xboxsize);
834 if (psr->yboxsize == 0)
836 psr->yboxsize = size/mat[0].nx;
837 printf("Set the y-size of the box to %.3f\n", psr->yboxsize);
841 for (i = 0; (i < nmat); i++)
843 if (mat[i].nmap > nmap1)
852 printf("Selected legend of matrix # %d for display\n", leg);
857 for (i = 0; (i < nmat); i++)
859 if (mat2[i].nmap > nmap2)
861 nmap2 = mat2[i].nmap;
868 printf("Selected legend of matrix # %d for second display\n", leg);
871 if ( (mat[0].legend[0] == 0) && psr->legend)
873 strcpy(mat[0].legend, psr->leglabel);
876 bTitle = bTitle && mat[nmat-1].title[0];
877 bTitleOnce = bTitleOnce && mat[nmat-1].title[0];
878 psr->bTitle = bTitle;
879 psr->bTitleOnce = bTitleOnce;
880 psr->bYonce = bYonce;
882 /* Set up size of box for nice colors */
883 box_dim(nmat, mat, mat2, psr, elegend, bFrame, &w, &h, &dw, &dh);
885 /* Set up bounding box */
899 out = ps_open(outf, 0, 0, x, y);
900 ps_linewidth(out, psr->linewidth);
901 ps_init_rgb_box(out, psr->xboxsize, psr->yboxsize);
902 ps_init_rgb_nbox(out, psr->xboxsize, psr->yboxsize);
903 ps_translate(out, psr->xoffs, psr->yoffs);
907 ps_comment(out, "Here starts the BOX drawing");
908 draw_boxes(out, x0, y0, w, nmat, mat, psr);
911 for (i = 0; (i < nmat); i++)
913 if (bTitle || (bTitleOnce && i == nmat-1) )
915 /* Print title, if any */
917 ps_strfont(out, psr->titfont, psr->titfontsize);
918 if (!mat2 || (strcmp(mat[i].title, mat2[i].title) == 0))
920 strcpy(buf, mat[i].title);
924 sprintf(buf, "%s / %s", mat[i].title, mat2[i].title);
926 ps_ctext(out, x0+w/2, y0+box_height(&(mat[i]), psr)+psr->titfontsize,
929 sprintf(buf, "Here starts the filling of box #%d", i);
930 ps_comment(out, buf);
931 for (x = 0; (x < mat[i].nx); x++)
936 xx = x0+x*psr->xboxsize;
937 ps_moveto(out, xx, y0);
939 bMap1 = (!mat2 || (x < y || (x == y && bFirstDiag)));
940 if ((bDiag) || (x != y))
942 col = mat[i].matrix[x][y];
948 for (nexty = 1; (nexty <= mat[i].ny); nexty++)
950 bNextMap1 = (!mat2 || (x < nexty || (x == nexty && bFirstDiag)));
951 /* TRUE: upper left -> map1 */
952 /* FALSE: lower right -> map2 */
953 if ((nexty == mat[i].ny) || (!bDiag && (x == nexty)))
959 nextcol = mat[i].matrix[x][nexty];
961 if ( (nexty == mat[i].ny) || (col != nextcol) || (bMap1 != bNextMap1) )
967 ps_rgb_nbox(out, &(mat[i].map[col].rgb), nexty-y);
971 ps_rgb_nbox(out, &(mat2[i].map[col].rgb), nexty-y);
976 ps_moverel(out, 0, psr->yboxsize);
984 y0 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
987 if (psr->X.lineatzero || psr->Y.lineatzero)
989 /* reset y0 for first box */
991 ps_comment(out, "Here starts the zero lines drawing");
992 draw_zerolines(out, x0, y0, w, nmat, mat, psr);
995 if (elegend != elNone)
997 ps_comment(out, "Now it's legend time!");
998 ps_linewidth(out, psr->linewidth);
999 if (mat2 == NULL || elegend != elSecond)
1001 bDiscrete = mat[0].bDiscrete;
1002 legend = mat[0].legend;
1008 bDiscrete = mat2[0].bDiscrete;
1009 legend = mat2[0].legend;
1015 leg_discrete(out, psr->legfontsize, DDD, legend,
1016 psr->legfontsize, psr->legfont, leg_nmap, leg_map);
1020 if (elegend != elBoth)
1022 leg_continuous(out, x0+w/2, w/2, DDD, legend,
1023 psr->legfontsize, psr->legfont, leg_nmap, leg_map,
1028 leg_bicontinuous(out, x0+w/2, w, DDD, mat[0].legend, mat2[0].legend,
1029 psr->legfontsize, psr->legfont, nmap1, map1, nmap2, map2);
1032 ps_comment(out, "Were there, dude");
1038 void make_axis_labels(int nmat, t_matrix *mat)
1042 for (i = 0; (i < nmat); i++)
1044 /* Make labels for x axis */
1045 if (mat[i].axis_x == NULL)
1047 snew(mat[i].axis_x, mat[i].nx);
1048 for (j = 0; (j < mat[i].nx); j++)
1050 mat[i].axis_x[j] = j;
1053 /* Make labels for y axis */
1054 if (mat[i].axis_y == NULL)
1056 snew(mat[i].axis_y, mat[i].ny);
1057 for (j = 0; (j < mat[i].ny); j++)
1059 mat[i].axis_y[j] = j;
1065 void prune_mat(int nmat, t_matrix *mat, t_matrix *mat2, int skip)
1067 int i, x, y, xs, ys;
1069 for (i = 0; i < nmat; i++)
1071 fprintf(stderr, "converting %dx%d matrix to %dx%d\n",
1072 mat[i].nx, mat[i].ny,
1073 (mat[i].nx+skip-1)/skip, (mat[i].ny+skip-1)/skip);
1074 /* walk through matrix */
1076 for (x = 0; (x < mat[i].nx); x++)
1080 mat[i].axis_x[xs] = mat[i].axis_x[x];
1083 mat2[i].axis_x[xs] = mat2[i].axis_x[x];
1086 for (y = 0; (y < mat[i].ny); y++)
1090 mat[i].axis_y[ys] = mat[i].axis_y[y];
1093 mat2[i].axis_y[ys] = mat2[i].axis_y[y];
1098 mat[i].matrix[xs][ys] = mat[i].matrix[x][y];
1101 mat2[i].matrix[xs][ys] = mat2[i].matrix[x][y];
1109 /* adjust parameters */
1110 mat[i].nx = (mat[i].nx+skip-1)/skip;
1111 mat[i].ny = (mat[i].ny+skip-1)/skip;
1114 mat2[i].nx = (mat2[i].nx+skip-1)/skip;
1115 mat2[i].ny = (mat2[i].ny+skip-1)/skip;
1120 void zero_lines(int nmat, t_matrix *mat, t_matrix *mat2)
1125 for (i = 0; i < nmat; i++)
1127 for (m = 0; m < (mat2 ? 2 : 1); m++)
1137 for (x = 0; x < mats[i].nx-1; x++)
1139 if (abs(mats[i].axis_x[x+1]) < 1e-5)
1141 for (y = 0; y < mats[i].ny; y++)
1143 mats[i].matrix[x][y] = 0;
1147 for (y = 0; y < mats[i].ny-1; y++)
1149 if (abs(mats[i].axis_y[y+1]) < 1e-5)
1151 for (x = 0; x < mats[i].nx; x++)
1153 mats[i].matrix[x][y] = 0;
1161 void write_combined_matrix(int ecombine, const char *fn,
1162 int nmat, t_matrix *mat1, t_matrix *mat2,
1163 real *cmin, real *cmax)
1165 int i, j, k, nlevels;
1166 t_mapping *map = NULL;
1168 real **rmat1, **rmat2;
1171 out = ffopen(fn, "w");
1172 for (k = 0; k < nmat; k++)
1174 if (mat2[k].nx != mat1[k].nx || mat2[k].ny != mat1[k].ny)
1176 gmx_fatal(FARGS, "Size of frame %d in 1st (%dx%d) and 2nd matrix (%dx%d) do"
1177 " not match.\n", k, mat1[k].nx, mat1[k].ny, mat2[k].nx, mat2[k].ny);
1179 printf("Combining two %dx%d matrices\n", mat1[k].nx, mat1[k].ny);
1180 rmat1 = matrix2real(&mat1[k], NULL);
1181 rmat2 = matrix2real(&mat2[k], NULL);
1182 if (NULL == rmat1 || NULL == rmat2)
1184 gmx_fatal(FARGS, "Could not extract real data from %s xpm matrices. Note that, e.g.,\n"
1185 "g_rms and g_mdmat provide such data, but not do_dssp.\n",
1186 (NULL == rmat1 && NULL == rmat2) ? "both" : "one of the" );
1190 for (j = 0; j < mat1[k].ny; j++)
1192 for (i = 0; i < mat1[k].nx; i++)
1196 case ecAdd: rmat1[i][j] += rmat2[i][j]; break;
1197 case ecSub: rmat1[i][j] -= rmat2[i][j]; break;
1198 case ecMult: rmat1[i][j] *= rmat2[i][j]; break;
1199 case ecDiv: rmat1[i][j] /= rmat2[i][j]; break;
1201 gmx_fatal(FARGS, "No such combination rule %d for matrices", ecombine);
1203 rlo = min(rlo, rmat1[i][j]);
1204 rhi = max(rhi, rmat1[i][j]);
1215 nlevels = max(mat1[k].nmap, mat2[k].nmap);
1219 "combination results in uniform matrix (%g), no output\n", rhi);
1222 else if (rlo>=0 || rhi<=0)
1223 write_xpm(out, mat1[k].flags, mat1[k].title, mat1[k].legend,
1224 mat1[k].label_x, mat1[k].label_y,
1225 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1226 rmat1, rlo, rhi, rhi<=0?red:white, rhi<=0?white:blue,
1229 write_xpm3(out, mat2[k].flags, mat1[k].title, mat1[k].legend,
1230 mat1[k].label_x, mat1[k].label_y,
1231 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1232 rmat1, rlo, 0, rhi, red, white, blue, &nlevels);
1236 write_xpm(out, mat1[k].flags, mat1[k].title, mat1[k].legend,
1237 mat1[k].label_x, mat1[k].label_y,
1238 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1239 rmat1, rlo, rhi, white, black, &nlevels);
1245 void do_mat(int nmat, t_matrix *mat, t_matrix *mat2,
1246 gmx_bool bFrame, gmx_bool bZeroLine, gmx_bool bDiag, gmx_bool bFirstDiag, gmx_bool bTitle,
1247 gmx_bool bTitleOnce, gmx_bool bYonce, int elegend,
1248 real size, real boxx, real boxy,
1249 const char *epsfile, const char *xpmfile, const char *m2p,
1250 const char *m2pout, int skip, int mapoffset)
1256 for (k = 0; (k < nmat); k++)
1258 if ((mat2[k].nx != mat[k].nx) || (mat2[k].ny != mat[k].ny))
1260 gmx_fatal(FARGS, "WAKE UP!! Size of frame %d in 2nd matrix file (%dx%d) does not match size of 1st matrix (%dx%d) or the other way around.\n",
1261 k, mat2[k].nx, mat2[k].ny, mat[k].nx, mat[k].ny);
1263 for (j = 0; (j < mat[k].ny); j++)
1265 for (i = bFirstDiag ? j+1 : j; (i < mat[k].nx); i++)
1267 mat[k].matrix[i][j] = mat2[k].matrix[i][j];
1272 for (i = 0; (i < nmat); i++)
1274 fprintf(stderr, "Matrix %d is %d x %d\n", i, mat[i].nx, mat[i].ny);
1277 make_axis_labels(nmat, mat);
1281 prune_mat(nmat, mat, mat2, skip);
1286 zero_lines(nmat, mat, mat);
1289 if (epsfile != NULL)
1291 ps_mat(epsfile, nmat, mat, mat2, bFrame, bDiag, bFirstDiag,
1292 bTitle, bTitleOnce, bYonce, elegend,
1293 size, boxx, boxy, m2p, m2pout, mapoffset);
1295 if (xpmfile != NULL)
1297 xpm_mat(xpmfile, nmat, mat, mat2, bDiag, bFirstDiag);
1301 void gradient_map(rvec grad, int nmap, t_mapping map[])
1306 for (i = 0; i < nmap; i++)
1309 map[i].rgb.r = 1-c*(1-grad[XX]);
1310 map[i].rgb.g = 1-c*(1-grad[YY]);
1311 map[i].rgb.b = 1-c*(1-grad[ZZ]);
1315 void gradient_mat(rvec grad, int nmat, t_matrix mat[])
1319 for (m = 0; m < nmat; m++)
1321 gradient_map(grad, mat[m].nmap, mat[m].map);
1325 void rainbow_map(gmx_bool bBlue, int nmap, t_mapping map[])
1330 for (i = 0; i < nmap; i++)
1332 c = (map[i].rgb.r + map[i].rgb.g + map[i].rgb.b)/3;
1341 if (c <= 0.25) /* 0-0.25 */
1344 g = pow(4*c, 0.666);
1347 else if (c <= 0.5) /* 0.25-0.5 */
1351 b = pow(2-4*c, 0.666);
1353 else if (c <= 0.75) /* 0.5-0.75 */
1355 r = pow(4*c-2, 0.666);
1362 g = pow(4-4*c, 0.666);
1371 void rainbow_mat(gmx_bool bBlue, int nmat, t_matrix mat[])
1375 for (m = 0; m < nmat; m++)
1377 rainbow_map(bBlue, mat[m].nmap, mat[m].map);
1381 int gmx_xpm2ps(int argc, char *argv[])
1383 const char *desc[] = {
1384 "[THISMODULE] makes a beautiful color plot of an XPixelMap file.",
1385 "Labels and axis can be displayed, when they are supplied",
1386 "in the correct matrix format.",
1387 "Matrix data may be generated by programs such as [gmx-do_dssp], [gmx-rms] or",
1388 "[gmx-mdmat].[PAR]",
1389 "Parameters are set in the [TT].m2p[tt] file optionally supplied with",
1390 "[TT]-di[tt]. Reasonable defaults are provided. Settings for the [IT]y[it]-axis",
1391 "default to those for the [IT]x[it]-axis. Font names have a defaulting hierarchy:",
1392 "titlefont -> legendfont; titlefont -> (xfont -> yfont -> ytickfont)",
1393 "-> xtickfont, e.g. setting titlefont sets all fonts, setting xfont",
1394 "sets yfont, ytickfont and xtickfont.[PAR]",
1395 "When no [TT].m2p[tt] file is supplied, many settings are taken from",
1396 "command line options. The most important option is [TT]-size[tt],",
1397 "which sets the size of the whole matrix in postscript units.",
1398 "This option can be overridden with the [TT]-bx[tt] and [TT]-by[tt]",
1399 "options (and the corresponding parameters in the [TT].m2p[tt] file),",
1400 "which set the size of a single matrix element.[PAR]",
1401 "With [TT]-f2[tt] a second matrix file can be supplied. Both matrix",
1402 "files will be read simultaneously and the upper left half of the",
1403 "first one ([TT]-f[tt]) is plotted together with the lower right",
1404 "half of the second one ([TT]-f2[tt]). The diagonal will contain",
1405 "values from the matrix file selected with [TT]-diag[tt].",
1406 "Plotting of the diagonal values can be suppressed altogether by",
1407 "setting [TT]-diag[tt] to [TT]none[tt].",
1408 "In this case, a new color map will be generated with",
1409 "a red gradient for negative numbers and a blue for positive.",
1410 "If the color coding and legend labels of both matrices are identical,",
1411 "only one legend will be displayed, else two separate legends are",
1413 "With [TT]-combine[tt], an alternative operation can be selected",
1414 "to combine the matrices. The output range is automatically set",
1415 "to the actual range of the combined matrix. This can be overridden",
1416 "with [TT]-cmin[tt] and [TT]-cmax[tt].[PAR]",
1417 "[TT]-title[tt] can be set to [TT]none[tt] to suppress the title, or to",
1418 "[TT]ylabel[tt] to show the title in the Y-label position (alongside",
1419 "the [IT]y[it]-axis).[PAR]",
1420 "With the [TT]-rainbow[tt] option, dull grayscale matrices can be turned",
1421 "into attractive color pictures.[PAR]",
1422 "Merged or rainbowed matrices can be written to an XPixelMap file with",
1423 "the [TT]-xpm[tt] option."
1427 const char *fn, *epsfile = NULL, *xpmfile = NULL;
1428 int i, nmat, nmat2, etitle, elegend, ediag, erainbow, ecombine;
1429 t_matrix *mat = NULL, *mat2 = NULL;
1430 gmx_bool bTitle, bTitleOnce, bDiag, bFirstDiag, bGrad;
1431 static gmx_bool bFrame = TRUE, bZeroLine = FALSE, bYonce = FALSE, bAdd = FALSE;
1432 static real size = 400, boxx = 0, boxy = 0, cmin = 0, cmax = 0;
1433 static rvec grad = {0, 0, 0};
1435 etSel, etTop, etOnce, etYlabel, etNone, etNR
1437 const char *title[] = { NULL, "top", "once", "ylabel", "none", NULL };
1438 /* MUST correspond to enum elXxx as defined at top of file */
1439 const char *legend[] = { NULL, "both", "first", "second", "none", NULL };
1441 edSel, edFirst, edSecond, edNone, edNR
1443 const char *diag[] = { NULL, "first", "second", "none", NULL };
1445 erSel, erNo, erBlue, erRed, erNR
1447 const char *rainbow[] = { NULL, "no", "blue", "red", NULL };
1448 /* MUST correspond to enum ecXxx as defined at top of file */
1449 const char *combine[] = {
1450 NULL, "halves", "add", "sub", "mult", "div", NULL
1452 static int skip = 1, mapoffset = 0;
1454 { "-frame", FALSE, etBOOL, {&bFrame},
1455 "Display frame, ticks, labels, title and legend" },
1456 { "-title", FALSE, etENUM, {title}, "Show title at" },
1457 { "-yonce", FALSE, etBOOL, {&bYonce}, "Show y-label only once" },
1458 { "-legend", FALSE, etENUM, {legend}, "Show legend" },
1459 { "-diag", FALSE, etENUM, {diag}, "Diagonal" },
1460 { "-size", FALSE, etREAL, {&size},
1461 "Horizontal size of the matrix in ps units" },
1462 { "-bx", FALSE, etREAL, {&boxx},
1463 "Element x-size, overrides [TT]-size[tt] (also y-size when [TT]-by[tt] is not set)" },
1464 { "-by", FALSE, etREAL, {&boxy}, "Element y-size" },
1465 { "-rainbow", FALSE, etENUM, {rainbow},
1466 "Rainbow colors, convert white to" },
1467 { "-gradient", FALSE, etRVEC, {grad},
1468 "Re-scale colormap to a smooth gradient from white {1,1,1} to {r,g,b}" },
1469 { "-skip", FALSE, etINT, {&skip},
1470 "only write out every nr-th row and column" },
1471 { "-zeroline", FALSE, etBOOL, {&bZeroLine},
1472 "insert line in [TT].xpm[tt] matrix where axis label is zero"},
1473 { "-legoffset", FALSE, etINT, {&mapoffset},
1474 "Skip first N colors from [TT].xpm[tt] file for the legend" },
1475 { "-combine", FALSE, etENUM, {combine}, "Combine two matrices" },
1476 { "-cmin", FALSE, etREAL, {&cmin}, "Minimum for combination output" },
1477 { "-cmax", FALSE, etREAL, {&cmax}, "Maximum for combination output" }
1479 #define NPA asize(pa)
1481 { efXPM, "-f", NULL, ffREAD },
1482 { efXPM, "-f2", "root2", ffOPTRD },
1483 { efM2P, "-di", NULL, ffLIBOPTRD },
1484 { efM2P, "-do", "out", ffOPTWR },
1485 { efEPS, "-o", NULL, ffOPTWR },
1486 { efXPM, "-xpm", NULL, ffOPTWR }
1488 #define NFILE asize(fnm)
1490 if (!parse_common_args(&argc, argv, PCA_CAN_VIEW,
1491 NFILE, fnm, NPA, pa,
1492 asize(desc), desc, 0, NULL, &oenv))
1497 etitle = nenum(title);
1498 elegend = nenum(legend);
1499 ediag = nenum(diag);
1500 erainbow = nenum(rainbow);
1501 ecombine = nenum(combine);
1502 bGrad = opt2parg_bSet("-gradient", NPA, pa);
1503 for (i = 0; i < DIM; i++)
1505 if (grad[i] < 0 || grad[i] > 1)
1507 gmx_fatal(FARGS, "RGB value %g out of range (0.0-1.0)", grad[i]);
1516 epsfile = ftp2fn_null(efEPS, NFILE, fnm);
1517 xpmfile = opt2fn_null("-xpm", NFILE, fnm);
1518 if (epsfile == NULL && xpmfile == NULL)
1520 if (ecombine != ecHalves)
1522 xpmfile = opt2fn("-xpm", NFILE, fnm);
1526 epsfile = ftp2fn(efEPS, NFILE, fnm);
1529 if (ecombine != ecHalves && epsfile)
1532 "WARNING: can only write result of arithmetic combination "
1533 "of two matrices to .xpm file\n"
1534 " file %s will not be written\n", epsfile);
1538 bDiag = ediag != edNone;
1539 bFirstDiag = ediag != edSecond;
1541 fn = opt2fn("-f", NFILE, fnm);
1542 nmat = read_xpm_matrix(fn, &mat);
1543 fprintf(stderr, "There %s %d matri%s in %s\n", (nmat > 1) ? "are" : "is", nmat, (nmat > 1) ? "ces" : "x", fn);
1544 fn = opt2fn_null("-f2", NFILE, fnm);
1547 nmat2 = read_xpm_matrix(fn, &mat2);
1548 fprintf(stderr, "There %s %d matri%s in %s\n", (nmat2 > 1) ? "are" : "is", nmat2, (nmat2 > 1) ? "ces" : "x", fn);
1551 fprintf(stderr, "Different number of matrices, using the smallest number.\n");
1552 nmat = nmat2 = min(nmat, nmat2);
1557 if (ecombine != ecHalves)
1560 "WARNING: arithmetic matrix combination selected (-combine), "
1561 "but no second matrix (-f2) supplied\n"
1562 " no matrix combination will be performed\n");
1567 bTitle = etitle == etTop;
1568 bTitleOnce = etitle == etOnce;
1569 if (etitle == etYlabel)
1571 for (i = 0; (i < nmat); i++)
1573 strcpy(mat[i].label_y, mat[i].title);
1576 strcpy(mat2[i].label_y, mat2[i].title);
1582 gradient_mat(grad, nmat, mat);
1585 gradient_mat(grad, nmat2, mat2);
1588 if (erainbow != erNo)
1590 rainbow_mat(erainbow == erBlue, nmat, mat);
1593 rainbow_mat(erainbow == erBlue, nmat2, mat2);
1597 if ((mat2 == NULL) && (elegend != elNone))
1602 if (ecombine && ecombine != ecHalves)
1604 write_combined_matrix(ecombine, xpmfile, nmat, mat, mat2,
1605 opt2parg_bSet("-cmin", NPA, pa) ? &cmin : NULL,
1606 opt2parg_bSet("-cmax", NPA, pa) ? &cmax : NULL);
1610 do_mat(nmat, mat, mat2, bFrame, bZeroLine, bDiag, bFirstDiag,
1611 bTitle, bTitleOnce, bYonce,
1612 elegend, size, boxx, boxy, epsfile, xpmfile,
1613 opt2fn_null("-di", NFILE, fnm), opt2fn_null("-do", NFILE, fnm), skip,
1617 view_all(oenv, NFILE, fnm);