?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Introduce SimulatorBuilder
[gromacs.git] / src / gromacs / gmxana / gmx_xpm2ps.cpp
blob8a881e1640d1266f848cfd9d34a76c02d6452b1d
1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, 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.
10 *
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.
15 *
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.
20 *
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.
25 *
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.
33 *
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.
36 */
37 #include "gmxpre.h"
38
39 #include <cassert>
40 #include <cmath>
41 #include <cstdio>
42 #include <cstdlib>
43 #include <cstring>
44
45 #include <algorithm>
46 #include <string>
47
48 #include "gromacs/commandline/pargs.h"
49 #include "gromacs/commandline/viewit.h"
50 #include "gromacs/fileio/matio.h"
51 #include "gromacs/fileio/readinp.h"
52 #include "gromacs/fileio/trxio.h"
53 #include "gromacs/fileio/warninp.h"
54 #include "gromacs/fileio/writeps.h"
55 #include "gromacs/gmxana/gmx_ana.h"
56 #include "gromacs/utility/arraysize.h"
57 #include "gromacs/utility/cstringutil.h"
58 #include "gromacs/utility/exceptions.h"
59 #include "gromacs/utility/fatalerror.h"
60 #include "gromacs/utility/filestream.h"
61 #include "gromacs/utility/futil.h"
62 #include "gromacs/utility/gmxassert.h"
63 #include "gromacs/utility/smalloc.h"
64 #include "gromacs/utility/stringutil.h"
65
66 #define FUDGE 1.2
67 #define DDD 2
68
69 typedef struct {
70 real major;
71 real minor;
72 real offset;
73 gmx_bool first;
74 int lineatzero;
75 real majorticklen;
76 real minorticklen;
77 char label[STRLEN];
78 real fontsize;
79 char font[STRLEN];
80 real tickfontsize;
81 char tickfont[STRLEN];
82 } t_axisdef;
83
84 typedef struct {
85 int bw;
86 real linewidth;
87 real xoffs, yoffs;
88 gmx_bool bTitle;
89 gmx_bool bTitleOnce;
90 gmx_bool bYonce;
91 real titfontsize;
92 char titfont[STRLEN];
93 gmx_bool legend;
94 real legfontsize;
95 char legfont[STRLEN];
96 char leglabel[STRLEN];
97 char leg2label[STRLEN];
98 real xboxsize;
99 real yboxsize;
100 real boxspacing;
101 real boxlinewidth;
102 real ticklinewidth;
103 real zerolinewidth;
104 t_axisdef X, Y;
105 } t_psrec;
106
107 /* MUST correspond to char *legend[] in main() */
108 enum {
109 elSel, elBoth, elFirst, elSecond, elNone, elNR
110 };
111
112 /* MUST correspond to char *combine[] in main() */
113 enum {
114 ecSel, ecHalves, ecAdd, ecSub, ecMult, ecDiv, ecNR
115 };
116
117 static void get_params(const char *mpin, const char *mpout, t_psrec *psr)
118 {
119 static const char *gmx_bools[BOOL_NR+1] = { "no", "yes", nullptr };
120 /* this must correspond to t_rgb *linecolors[] below */
121 static const char *colors[] = { "none", "black", "white", nullptr };
122 warninp_t wi;
123 std::vector<t_inpfile> inp;
124
125 wi = init_warning(FALSE, 0);
126
127 if (mpin != nullptr)
128 {
129 std::string libmpin = gmx::findLibraryFile(mpin);
130 gmx::TextInputFile stream(libmpin);
131 inp = read_inpfile(&stream, libmpin.c_str(), wi);
132 }
133 else
134 {
135 inp.clear();
136 }
137
138 psr->bw = get_eenum(&inp, "black&white", gmx_bools);
139 psr->linewidth = get_ereal(&inp, "linewidth", 1.0, wi);
140 setStringEntry(&inp, "titlefont", psr->titfont, "Helvetica");
141 psr->titfontsize = get_ereal(&inp, "titlefontsize", 20.0, wi);
142 psr->legend = (get_eenum(&inp, "legend", gmx_bools) != 0);
143 setStringEntry(&inp, "legendfont", psr->legfont, psr->titfont);
144 setStringEntry(&inp, "legendlabel", psr->leglabel, "");
145 setStringEntry(&inp, "legend2label", psr->leg2label, psr->leglabel);
146 psr->legfontsize = get_ereal(&inp, "legendfontsize", 14.0, wi);
147 psr->xboxsize = get_ereal(&inp, "xbox", 0.0, wi);
148 psr->yboxsize = get_ereal(&inp, "ybox", 0.0, wi);
149 psr->boxspacing = get_ereal(&inp, "matrixspacing", 20.0, wi);
150 psr->xoffs = get_ereal(&inp, "xoffset", 0.0, wi);
151 psr->yoffs = get_ereal(&inp, "yoffset", psr->xoffs, wi);
152 psr->boxlinewidth = get_ereal(&inp, "boxlinewidth", psr->linewidth, wi);
153 psr->ticklinewidth = get_ereal(&inp, "ticklinewidth", psr->linewidth, wi);
154 psr->zerolinewidth = get_ereal(&inp, "zerolinewidth", psr->ticklinewidth, wi);
155 psr->X.lineatzero = get_eenum(&inp, "x-lineat0value", colors);
156 psr->X.major = get_ereal(&inp, "x-major", 1, wi);
157 psr->X.minor = get_ereal(&inp, "x-minor", 1, wi);
158 psr->X.offset = get_ereal(&inp, "x-firstmajor", 0.0, wi);
159 psr->X.first = (get_eenum(&inp, "x-majorat0", gmx_bools) != 0);
160 psr->X.majorticklen = get_ereal(&inp, "x-majorticklen", 8.0, wi);
161 psr->X.minorticklen = get_ereal(&inp, "x-minorticklen", 4.0, wi);
162 setStringEntry(&inp, "x-label", psr->X.label, "");
163 psr->X.fontsize = get_ereal(&inp, "x-fontsize", 16.0, wi);
164 setStringEntry(&inp, "x-font", psr->X.font, psr->titfont);
165 psr->X.tickfontsize = get_ereal(&inp, "x-tickfontsize", 10.0, wi);
166 setStringEntry(&inp, "x-tickfont", psr->X.tickfont, psr->X.font);
167 psr->Y.lineatzero = get_eenum(&inp, "y-lineat0value", colors);
168 psr->Y.major = get_ereal(&inp, "y-major", psr->X.major, wi);
169 psr->Y.minor = get_ereal(&inp, "y-minor", psr->X.minor, wi);
170 psr->Y.offset = get_ereal(&inp, "y-firstmajor", psr->X.offset, wi);
171 psr->Y.first = (get_eenum(&inp, "y-majorat0", gmx_bools) != 0);
172 psr->Y.majorticklen = get_ereal(&inp, "y-majorticklen", psr->X.majorticklen, wi);
173 psr->Y.minorticklen = get_ereal(&inp, "y-minorticklen", psr->X.minorticklen, wi);
174 setStringEntry(&inp, "y-label", psr->Y.label, psr->X.label);
175 psr->Y.fontsize = get_ereal(&inp, "y-fontsize", psr->X.fontsize, wi);
176 setStringEntry(&inp, "y-font", psr->Y.font, psr->X.font);
177 psr->Y.tickfontsize = get_ereal(&inp, "y-tickfontsize", psr->X.tickfontsize, wi);
178 setStringEntry(&inp, "y-tickfont", psr->Y.tickfont, psr->Y.font);
179
180 check_warning_error(wi, FARGS);
181
182 if (mpout != nullptr)
183 {
184 gmx::TextOutputFile stream(mpout);
185 write_inpfile(&stream, mpout, &inp, TRUE, WriteMdpHeader::yes, wi);
186 stream.close();
187 }
188
189 done_warning(wi, FARGS);
190 }
191
192 static t_rgb black = { 0, 0, 0 };
193 static t_rgb white = { 1, 1, 1 };
194 #define BLACK (&black)
195 /* this must correspond to *colors[] in get_params */
196 static t_rgb *linecolors[] = { nullptr, &black, &white, nullptr };
197
198 static gmx_bool diff_maps(int nmap1, t_mapping *map1, int nmap2, t_mapping *map2)
199 {
200 int i;
201 gmx_bool bDiff, bColDiff = FALSE;
202
203 if (nmap1 != nmap2)
204 {
205 bDiff = TRUE;
206 }
207 else
208 {
209 bDiff = FALSE;
210 for (i = 0; i < nmap1; i++)
211 {
212 if (!matelmt_cmp(map1[i].code, map2[i].code))
213 {
214 bDiff = TRUE;
215 }
216 if (std::strcmp(map1[i].desc, map2[i].desc) != 0)
217 {
218 bDiff = TRUE;
219 }
220 if ((map1[i].rgb.r != map2[i].rgb.r) ||
221 (map1[i].rgb.g != map2[i].rgb.g) ||
222 (map1[i].rgb.b != map2[i].rgb.b))
223 {
224 bColDiff = TRUE;
225 }
226 }
227 if (!bDiff && bColDiff)
228 {
229 fprintf(stderr, "Warning: two colormaps differ only in RGB value, using one colormap.\n");
230 }
231 }
232
233 return bDiff;
234 }
235
236 static void leg_discrete(t_psdata ps, real x0, real y0, char *label,
237 real fontsize, char *font, int nmap, t_mapping map[])
238 {
239 int i;
240 real yhh;
241 real boxhh;
242
243 boxhh = fontsize+DDD;
244 /* LANDSCAPE */
245 ps_rgb(ps, BLACK);
246 ps_strfont(ps, font, fontsize);
247 yhh = y0+fontsize+3*DDD;
248 if (std::strlen(label) > 0)
249 {
250 ps_ctext(ps, x0, yhh, label, eXLeft);
251 }
252 ps_moveto(ps, x0, y0);
253 for (i = 0; (i < nmap); i++)
254 {
255 ps_setorigin(ps);
256 ps_rgb(ps, &(map[i].rgb));
257 ps_fillbox(ps, DDD, DDD, DDD+fontsize, boxhh-DDD);
258 ps_rgb(ps, BLACK);
259 ps_box(ps, DDD, DDD, DDD+fontsize, boxhh-DDD);
260 ps_ctext(ps, boxhh+2*DDD, fontsize/3, map[i].desc, eXLeft);
261 ps_unsetorigin(ps);
262 ps_moverel(ps, DDD, -fontsize/3);
263 }
264 }
265
266 static void leg_continuous(t_psdata ps, real x0, real x, real y0, char *label,
267 real fontsize, char *font,
268 int nmap, t_mapping map[],
269 int mapoffset)
270 {
271 int i;
272 real xx0;
273 real yhh, boxxh, boxyh;
274
275 boxyh = fontsize;
276 if (x < 8*fontsize)
277 {
278 x = 8*fontsize;
279 }
280 boxxh = x/(nmap-mapoffset);
281 if (boxxh > fontsize)
282 {
283 boxxh = fontsize;
284 }
285
286 GMX_RELEASE_ASSERT(map != nullptr, "NULL map array provided to leg_continuous()");
287
288 /* LANDSCAPE */
289 xx0 = x0-((nmap-mapoffset)*boxxh)/2.0;
290
291 for (i = 0; (i < nmap-mapoffset); i++)
292 {
293 ps_rgb(ps, &(map[i+mapoffset].rgb));
294 ps_fillbox(ps, xx0+i*boxxh, y0, xx0+(i+1)*boxxh, y0+boxyh);
295 }
296 ps_strfont(ps, font, fontsize);
297 ps_rgb(ps, BLACK);
298 ps_box(ps, xx0, y0, xx0+(nmap-mapoffset)*boxxh, y0+boxyh);
299
300 yhh = y0+boxyh+3*DDD;
301 ps_ctext(ps, xx0+boxxh/2, yhh, map[0].desc, eXCenter);
302 if (std::strlen(label) > 0)
303 {
304 ps_ctext(ps, x0, yhh, label, eXCenter);
305 }
306 ps_ctext(ps, xx0+((nmap-mapoffset)*boxxh)
307 - boxxh/2, yhh, map[nmap-1].desc, eXCenter);
308 }
309
310 static void leg_bicontinuous(t_psdata ps, real x0, real x, real y0, char *label1,
311 char *label2, real fontsize, char *font,
312 int nmap1, t_mapping map1[], int nmap2, t_mapping map2[])
313 {
314 real xx1, xx2, x1, x2;
315
316 x1 = x/(nmap1+nmap2)*nmap1; /* width of legend 1 */
317 x2 = x/(nmap1+nmap2)*nmap2; /* width of legend 2 */
318 xx1 = x0-(x2/2.0)-fontsize; /* center of legend 1 */
319 xx2 = x0+(x1/2.0)+fontsize; /* center of legend 2 */
320 x1 -= fontsize/2; /* adjust width */
321 x2 -= fontsize/2; /* adjust width */
322 leg_continuous(ps, xx1, x1, y0, label1, fontsize, font, nmap1, map1, 0);
323 leg_continuous(ps, xx2, x2, y0, label2, fontsize, font, nmap2, map2, 0);
324 }
325
326 static real box_height(t_matrix *mat, t_psrec *psr)
327 {
328 return mat->ny*psr->yboxsize;
329 }
330
331 static real box_dh(t_psrec *psr)
332 {
333 return psr->boxspacing;
334 }
335
336 #define IS_ONCE (i == nmat-1)
337 static real box_dh_top(gmx_bool bOnce, t_psrec *psr)
338 {
339 real dh;
340
341 if (psr->bTitle || (psr->bTitleOnce && bOnce) )
342 {
343 dh = 2*psr->titfontsize;
344 }
345 else
346 {
347 dh = 0;
348 }
349
350 return dh;
351 }
352
353 static gmx_bool box_do_all_x_maj_ticks(t_psrec *psr)
354 {
355 return (psr->boxspacing > (1.5*psr->X.majorticklen));
356 }
357
358 static gmx_bool box_do_all_x_min_ticks(t_psrec *psr)
359 {
360 return (psr->boxspacing > (1.5*psr->X.minorticklen));
361 }
362
363 static void draw_boxes(t_psdata ps, real x0, real y0, real w,
364 int nmat, t_matrix mat[], t_psrec *psr)
365 {
366 char buf[128];
367 char *mylab;
368 real xxx;
369 char **xtick, **ytick;
370 real xx, yy, dy, xx00, yy00, offset_x, offset_y;
371 int i, j, x, y, ntx, nty;
372 size_t strlength;
373
374 /* Only necessary when there will be no y-labels */
375 strlength = 0;
376
377 /* Draw the box */
378 ps_rgb(ps, BLACK);
379 ps_linewidth(ps, static_cast<int>(psr->boxlinewidth));
380 yy00 = y0;
381 for (i = 0; (i < nmat); i++)
382 {
383 dy = box_height(&(mat[i]), psr);
384 ps_box(ps, x0-1, yy00-1, x0+w+1, yy00+dy+1);
385 yy00 += dy+box_dh(psr)+box_dh_top(IS_ONCE, psr);
386 }
387
388 /* Draw the ticks on the axes */
389 ps_linewidth(ps, static_cast<int>(psr->ticklinewidth));
390 xx00 = x0-1;
391 yy00 = y0-1;
392 for (i = 0; (i < nmat); i++)
393 {
394 if (mat[i].flags & MAT_SPATIAL_X)
395 {
396 ntx = mat[i].nx + 1;
397 offset_x = 0.1;
398 }
399 else
400 {
401 ntx = mat[i].nx;
402 offset_x = 0.6;
403 }
404 if (mat[i].flags & MAT_SPATIAL_Y)
405 {
406 nty = mat[i].ny + 1;
407 offset_y = 0.1;
408 }
409 else
410 {
411 nty = mat[i].ny;
412 offset_y = 0.6;
413 }
414 snew(xtick, ntx);
415 for (j = 0; (j < ntx); j++)
416 {
417 sprintf(buf, "%g", mat[i].axis_x[j]);
418 xtick[j] = gmx_strdup(buf);
419 }
420 ps_strfont(ps, psr->X.tickfont, psr->X.tickfontsize);
421 for (x = 0; (x < ntx); x++)
422 {
423 xx = xx00 + (x + offset_x)*psr->xboxsize;
424 if ( ( bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.major) ||
425 (psr->X.first && (x == 0))) &&
426 ( (i == 0) || box_do_all_x_maj_ticks(psr) ) )
427 {
428 /* Longer tick marks */
429 ps_line (ps, xx, yy00, xx, yy00-psr->X.majorticklen);
430 /* Plot label on lowest graph only */
431 if (i == 0)
432 {
433 ps_ctext(ps, xx,
434 yy00-DDD-psr->X.majorticklen-psr->X.tickfontsize*0.8,
435 xtick[x], eXCenter);
436 }
437 }
438 else if (bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.minor) &&
439 ( (i == 0) || box_do_all_x_min_ticks(psr) ) )
440 {
441 /* Shorter tick marks */
442 ps_line(ps, xx, yy00, xx, yy00-psr->X.minorticklen);
443 }
444 else if (bRmod(mat[i].axis_x[x], psr->X.offset, psr->X.major) )
445 {
446 /* Even shorter marks, only each X.major */
447 ps_line(ps, xx, yy00, xx, yy00-(psr->boxspacing/2));
448 }
449 }
450 ps_strfont(ps, psr->Y.tickfont, psr->Y.tickfontsize);
451 snew(ytick, nty);
452 for (j = 0; (j < nty); j++)
453 {
454 sprintf(buf, "%g", mat[i].axis_y[j]);
455 ytick[j] = gmx_strdup(buf);
456 }
457
458 for (y = 0; (y < nty); y++)
459 {
460 yy = yy00 + (y + offset_y)*psr->yboxsize;
461 if (bRmod(mat[i].axis_y[y], psr->Y.offset, psr->Y.major) ||
462 (psr->Y.first && (y == 0)))
463 {
464 /* Major ticks */
465 strlength = std::max(strlength, std::strlen(ytick[y]));
466 ps_line (ps, xx00, yy, xx00-psr->Y.majorticklen, yy);
467 ps_ctext(ps, xx00-psr->Y.majorticklen-DDD,
468 yy-psr->Y.tickfontsize/3.0, ytick[y], eXRight);
469 }
470 else if (bRmod(mat[i].axis_y[y], psr->Y.offset, psr->Y.minor) )
471 {
472 /* Minor ticks */
473 ps_line(ps, xx00, yy, xx00-psr->Y.minorticklen, yy);
474 }
475 }
476 sfree(xtick);
477 sfree(ytick);
478
479 /* Label on Y-axis */
480 if (!psr->bYonce || i == nmat/2)
481 {
482 if (strlen(psr->Y.label) > 0)
483 {
484 mylab = psr->Y.label;
485 }
486 else
487 {
488 mylab = mat[i].label_y;
489 }
490 if (strlen(mylab) > 0)
491 {
492 ps_strfont(ps, psr->Y.font, psr->Y.fontsize);
493 ps_flip(ps, TRUE);
494 xxx = x0-psr->X.majorticklen-psr->X.tickfontsize*strlength-DDD;
495 ps_ctext(ps, yy00+box_height(&mat[i], psr)/2.0, 612.5-xxx,
496 mylab, eXCenter);
497 ps_flip(ps, FALSE);
498 }
499 }
500
501 yy00 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
502 }
503 /* Label on X-axis */
504 if (strlen(psr->X.label) > 0)
505 {
506 mylab = psr->X.label;
507 }
508 else
509 {
510 mylab = mat[0].label_x;
511 }
512 if (strlen(mylab) > 0)
513 {
514 ps_strfont(ps, psr->X.font, psr->X.fontsize);
515 ps_ctext(ps, x0+w/2, y0-DDD-psr->X.majorticklen-psr->X.tickfontsize*FUDGE-
516 psr->X.fontsize, mylab, eXCenter);
517 }
518 }
519
520 static void draw_zerolines(t_psdata out, real x0, real y0, real w,
521 int nmat, t_matrix mat[], t_psrec *psr)
522 {
523 real xx, yy, dy, xx00, yy00;
524 int i, x, y;
525
526 xx00 = x0-1.5;
527 yy00 = y0-1.5;
528 ps_linewidth(out, static_cast<int>(psr->zerolinewidth));
529 for (i = 0; (i < nmat); i++)
530 {
531 dy = box_height(&(mat[i]), psr);
532 /* mat[i].axis_x and _y were already set by draw_boxes */
533 if (psr->X.lineatzero)
534 {
535 ps_rgb(out, linecolors[psr->X.lineatzero]);
536 for (x = 0; (x < mat[i].nx); x++)
537 {
538 xx = xx00+(x+0.7)*psr->xboxsize;
539 /* draw lines whenever tick label almost zero (e.g. next trajectory) */
540 if (x != 0 && x < mat[i].nx-1 &&
541 std::abs(mat[i].axis_x[x]) <
542 0.1*std::abs(mat[i].axis_x[x+1]-mat[i].axis_x[x]) )
543 {
544 ps_line (out, xx, yy00, xx, yy00+dy+2);
545 }
546 }
547 }
548 if (psr->Y.lineatzero)
549 {
550 ps_rgb(out, linecolors[psr->Y.lineatzero]);
551 for (y = 0; (y < mat[i].ny); y++)
552 {
553 yy = yy00+(y+0.7)*psr->yboxsize;
554 /* draw lines whenever tick label almost zero (e.g. next trajectory) */
555 if (y != 0 && y < mat[i].ny-1 &&
556 std::abs(mat[i].axis_y[y]) <
557 0.1*std::abs(mat[i].axis_y[y+1]-mat[i].axis_y[y]) )
558 {
559 ps_line (out, xx00, yy, xx00+w+2, yy);
560 }
561 }
562 }
563 yy00 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
564 }
565 }
566
567 static void box_dim(int nmat, t_matrix mat[], t_matrix *mat2, t_psrec *psr,
568 int elegend, gmx_bool bFrame,
569 real *w, real *h, real *dw, real *dh)
570 {
571 int i, maxytick;
572 real ww, hh, dww, dhh;
573
574 hh = dww = dhh = 0;
575 maxytick = 0;
576
577 ww = 0;
578 for (i = 0; (i < nmat); i++)
579 {
580 ww = std::max(ww, mat[i].nx*psr->xboxsize);
581 hh += box_height(&(mat[i]), psr);
582 maxytick = std::max(maxytick, mat[i].nx);
583 }
584 if (bFrame)
585 {
586 if (mat[0].label_y[0])
587 {
588 dww += 2.0*(psr->Y.fontsize+DDD);
589 }
590 if (psr->Y.major > 0)
591 {
592 dww += psr->Y.majorticklen + DDD +
593 psr->Y.tickfontsize*(std::log(static_cast<real>(maxytick))/std::log(10.0));
594 }
595 else if (psr->Y.minor > 0)
596 {
597 dww += psr->Y.minorticklen;
598 }
599
600 if (mat[0].label_x[0])
601 {
602 dhh += psr->X.fontsize+2*DDD;
603 }
604 if ( /* fool emacs auto-indent */
605 (elegend == elBoth && (mat[0].legend[0] || (mat2 && mat2[0].legend[0]))) ||
606 (elegend == elFirst && mat[0].legend[0]) ||
607 (elegend == elSecond && (mat2 && mat2[0].legend[0])) )
608 {
609 dhh += 2*(psr->legfontsize*FUDGE+2*DDD);
610 }
611 else
612 {
613 dhh += psr->legfontsize*FUDGE+2*DDD;
614 }
615 if (psr->X.major > 0)
616 {
617 dhh += psr->X.tickfontsize*FUDGE+2*DDD+psr->X.majorticklen;
618 }
619 else if (psr->X.minor > 0)
620 {
621 dhh += psr->X.minorticklen;
622 }
623
624 hh += (nmat-1)*box_dh(psr);
625 hh += box_dh_top(TRUE, psr);
626 if (nmat > 1)
627 {
628 hh += (nmat-1)*box_dh_top(FALSE, psr);
629 }
630 }
631 *w = ww;
632 *h = hh;
633 *dw = dww;
634 *dh = dhh;
635 }
636
637 static int add_maps(t_mapping **newmap,
638 int nmap1, t_mapping map1[], int nmap2, t_mapping map2[])
639 {
640 static char mapper[] = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789!@#$%^&*()-_=+{}|;:',<.>/?";
641 int nsymbols;
642 int nmap, j, k;
643 t_mapping *map;
644
645 nsymbols = std::strlen(mapper);
646 nmap = nmap1+nmap2;
647 if (nmap > nsymbols*nsymbols)
648 {
649 gmx_fatal(FARGS, "Not enough symbols to merge the two colormaps\n");
650 }
651 printf("Combining colormaps of %d and %d elements into one of %d elements\n",
652 nmap1, nmap2, nmap);
653 snew(map, nmap);
654 for (j = 0; j < nmap1; j++)
655 {
656 map[j].code.c1 = mapper[j % nsymbols];
657 if (nmap > nsymbols)
658 {
659 map[j].code.c2 = mapper[j/nsymbols];
660 }
661 map[j].rgb.r = map1[j].rgb.r;
662 map[j].rgb.g = map1[j].rgb.g;
663 map[j].rgb.b = map1[j].rgb.b;
664 map[j].desc = map1[j].desc;
665 }
666 for (j = 0; j < nmap2; j++)
667 {
668 k = j+nmap1;
669 map[k].code.c1 = mapper[k % nsymbols];
670 if (nmap > nsymbols)
671 {
672 map[k].code.c2 = mapper[k/nsymbols];
673 }
674 map[k].rgb.r = map2[j].rgb.r;
675 map[k].rgb.g = map2[j].rgb.g;
676 map[k].rgb.b = map2[j].rgb.b;
677 map[k].desc = map2[j].desc;
678 }
679
680 *newmap = map;
681 return nmap;
682 }
683
684 static void xpm_mat(const char *outf, int nmat, t_matrix *mat, t_matrix *mat2,
685 gmx_bool bDiag, gmx_bool bFirstDiag)
686 {
687 FILE *out;
688 int i, x, y, col;
689 int nmap;
690 t_mapping *map = nullptr;
691
692 out = gmx_ffopen(outf, "w");
693
694 for (i = 0; i < nmat; i++)
695 {
696 if (!mat2 || !diff_maps(mat[i].nmap, mat[i].map, mat2[i].nmap, mat2[i].map))
697 {
698 write_xpm_m(out, mat[0]);
699 }
700 else
701 {
702 nmap = add_maps(&map, mat[i].nmap, mat[i].map, mat2[i].nmap, mat2[i].map);
703 for (x = 0; (x < mat[i].nx); x++)
704 {
705 for (y = 0; (y < mat[i].nx); y++)
706 {
707 if ((x < y) || ((x == y) && bFirstDiag)) /* upper left -> map1 */
708 {
709 col = mat[i].matrix[x][y];
710 }
711 else /* lower right -> map2 */
712 {
713 col = mat[i].nmap+mat[i].matrix[x][y];
714 }
715 if ((bDiag) || (x != y))
716 {
717 mat[i].matrix[x][y] = col;
718 }
719 else
720 {
721 mat[i].matrix[x][y] = 0;
722 }
723 }
724 }
725 sfree(mat[i].map);
726 mat[i].nmap = nmap;
727 mat[i].map = map;
728 if (std::strcmp(mat[i].title, mat2[i].title) != 0)
729 {
730 sprintf(mat[i].title+strlen(mat[i].title), " / %s", mat2[i].title);
731 }
732 if (std::strcmp(mat[i].legend, mat2[i].legend) != 0)
733 {
734 sprintf(mat[i].legend+strlen(mat[i].legend), " / %s", mat2[i].legend);
735 }
736 write_xpm_m(out, mat[i]);
737 }
738 }
739 gmx_ffclose(out);
740 }
741
742 static void tick_spacing(int n, real axis[], real offset, char axisnm,
743 real *major, real *minor)
744 {
745 real space;
746 gmx_bool bTryAgain;
747 int i, j, t, f = 0, ten;
748 #define NFACT 4
749 real major_fact[NFACT] = {5, 4, 2, 1};
750 real minor_fact[NFACT] = {5, 4, 4, 5};
751
752 /* start with interval between 10 matrix points: */
753 space = std::max(10*axis[1]-axis[0], axis[std::min(10, n-1)]-axis[0]);
754 /* get power of 10 */
755 ten = static_cast<int>(std::ceil(std::log(space)/std::log(10.0))-1);
756 bTryAgain = TRUE;
757 for (t = ten+2; t > ten-3 && bTryAgain; t--)
758 {
759 for (f = 0; f < NFACT && bTryAgain; f++)
760 {
761 space = std::pow(10.0_real, static_cast<real>(t)) * major_fact[f];
762 /* count how many ticks we would get: */
763 i = 0;
764 for (j = 0; j < n; j++)
765 {
766 if (bRmod(axis[j], offset, space) )
767 {
768 i++;
769 }
770 }
771 /* do we have a reasonable number of ticks ? */
772 bTryAgain = (i > std::min(10, n-1)) || (i < 5);
773 }
774 }
775 if (bTryAgain)
776 {
777 space = std::max(10*axis[1]-axis[0], axis[std::min(10, n-1)]-axis[0]);
778 fprintf(stderr, "Auto tick spacing failed for %c-axis, guessing %g\n",
779 axisnm, space);
780 }
781 *major = space;
782 *minor = space / minor_fact[(f > 0) ? f-1 : 0];
783 fprintf(stderr, "Auto tick spacing for %c-axis: major %g, minor %g\n",
784 axisnm, *major, *minor);
785 }
786
787 static void ps_mat(const char *outf, int nmat, t_matrix mat[], t_matrix mat2[],
788 gmx_bool bFrame, gmx_bool bDiag, gmx_bool bFirstDiag,
789 gmx_bool bTitle, gmx_bool bTitleOnce, gmx_bool bYonce, int elegend,
790 real size, real boxx, real boxy, const char *m2p, const char *m2pout,
791 int mapoffset)
792 {
793 char *legend;
794 t_psdata out;
795 t_psrec psrec, *psr;
796 int W, H;
797 int i, x, y, col, leg = 0;
798 real x0, y0, xx;
799 real w, h, dw, dh;
800 int nmap1 = 0, nmap2 = 0, leg_nmap;
801 t_mapping *map1 = nullptr, *map2 = nullptr, *leg_map;
802 gmx_bool bMap1, bNextMap1, bDiscrete;
803
804 try
805 {
806 get_params(m2p, m2pout, &psrec);
807 }
808 GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR;
809
810 psr = &psrec;
811
812 if (psr->X.major <= 0)
813 {
814 tick_spacing((mat[0].flags & MAT_SPATIAL_X) ? mat[0].nx + 1 : mat[0].nx,
815 mat[0].axis_x, psr->X.offset, 'X',
816 &(psr->X.major), &(psr->X.minor) );
817 }
818 if (psr->X.minor <= 0)
819 {
820 psr->X.minor = psr->X.major / 2;
821 }
822 if (psr->Y.major <= 0)
823 {
824 tick_spacing((mat[0].flags & MAT_SPATIAL_Y) ? mat[0].ny + 1 : mat[0].ny,
825 mat[0].axis_y, psr->Y.offset, 'Y',
826 &(psr->Y.major), &(psr->Y.minor) );
827 }
828 if (psr->Y.minor <= 0)
829 {
830 psr->Y.minor = psr->Y.major / 2;
831 }
832
833 if (boxx > 0)
834 {
835 psr->xboxsize = boxx;
836 psr->yboxsize = boxx;
837 }
838 if (boxy > 0)
839 {
840 psr->yboxsize = boxy;
841 }
842
843 if (psr->xboxsize == 0)
844 {
845 psr->xboxsize = size/mat[0].nx;
846 printf("Set the x-size of the box to %.3f\n", psr->xboxsize);
847 }
848 if (psr->yboxsize == 0)
849 {
850 psr->yboxsize = size/mat[0].nx;
851 printf("Set the y-size of the box to %.3f\n", psr->yboxsize);
852 }
853
854 nmap1 = 0;
855 for (i = 0; (i < nmat); i++)
856 {
857 if (mat[i].nmap > nmap1)
858 {
859 nmap1 = mat[i].nmap;
860 map1 = mat[i].map;
861 leg = i+1;
862 }
863 }
864 if (leg != 1)
865 {
866 printf("Selected legend of matrix # %d for display\n", leg);
867 }
868 if (mat2)
869 {
870 nmap2 = 0;
871 for (i = 0; (i < nmat); i++)
872 {
873 if (mat2[i].nmap > nmap2)
874 {
875 nmap2 = mat2[i].nmap;
876 map2 = mat2[i].map;
877 leg = i+1;
878 }
879 }
880 if (leg != 1)
881 {
882 printf("Selected legend of matrix # %d for second display\n", leg);
883 }
884 }
885 if ( (mat[0].legend[0] == 0) && psr->legend)
886 {
887 std::strcpy(mat[0].legend, psr->leglabel);
888 }
889
890 bTitle = bTitle && (mat[nmat-1].title[0] != 0);
891 bTitleOnce = bTitleOnce && (mat[nmat-1].title[0] != 0);
892 psr->bTitle = bTitle;
893 psr->bTitleOnce = bTitleOnce;
894 psr->bYonce = bYonce;
895
896 /* Set up size of box for nice colors */
897 box_dim(nmat, mat, mat2, psr, elegend, bFrame, &w, &h, &dw, &dh);
898
899 /* Set up bounding box */
900 W = static_cast<int>(w+dw);
901 H = static_cast<int>(h+dh);
902
903 /* Start box at */
904 x0 = dw;
905 y0 = dh;
906 x = static_cast<int>(W+psr->xoffs);
907 y = static_cast<int>(H+psr->yoffs);
908 if (bFrame)
909 {
910 x += 5*DDD;
911 y += 4*DDD;
912 }
913 out = ps_open(outf, 0, 0, x, y);
914 ps_linewidth(out, static_cast<int>(psr->linewidth));
915 ps_init_rgb_box(out, psr->xboxsize, psr->yboxsize);
916 ps_init_rgb_nbox(out, psr->xboxsize, psr->yboxsize);
917 ps_translate(out, psr->xoffs, psr->yoffs);
918
919 if (bFrame)
920 {
921 ps_comment(out, "Here starts the BOX drawing");
922 draw_boxes(out, x0, y0, w, nmat, mat, psr);
923 }
924
925 for (i = 0; (i < nmat); i++)
926 {
927 if (bTitle || (bTitleOnce && i == nmat-1) )
928 {
929 /* Print title, if any */
930 ps_rgb(out, BLACK);
931 ps_strfont(out, psr->titfont, psr->titfontsize);
932 std::string buf;
933 if (!mat2 || (std::strcmp(mat[i].title, mat2[i].title) == 0))
934 {
935 buf = mat[i].title;
936 }
937 else
938 {
939 buf = gmx::formatString("%s / %s", mat[i].title, mat2[i].title);
940 }
941 ps_ctext(out, x0+w/2, y0+box_height(&(mat[i]), psr)+psr->titfontsize,
942 buf.c_str(), eXCenter);
943 }
944 ps_comment(out, gmx::formatString("Here starts the filling of box #%d", i).c_str());
945 for (x = 0; (x < mat[i].nx); x++)
946 {
947 int nexty;
948 int nextcol;
949
950 xx = x0+x*psr->xboxsize;
951 ps_moveto(out, xx, y0);
952 y = 0;
953 bMap1 = ((mat2 == nullptr) || (x < y || (x == y && bFirstDiag)));
954 if ((bDiag) || (x != y))
955 {
956 col = mat[i].matrix[x][y];
957 }
958 else
959 {
960 col = -1;
961 }
962 for (nexty = 1; (nexty <= mat[i].ny); nexty++)
963 {
964 bNextMap1 = ((mat2 == nullptr) || (x < nexty || (x == nexty && bFirstDiag)));
965 /* TRUE: upper left -> map1 */
966 /* FALSE: lower right -> map2 */
967 if ((nexty == mat[i].ny) || (!bDiag && (x == nexty)))
968 {
969 nextcol = -1;
970 }
971 else
972 {
973 nextcol = mat[i].matrix[x][nexty];
974 }
975 if ( (nexty == mat[i].ny) || (col != nextcol) || (bMap1 != bNextMap1) )
976 {
977 if (col >= 0)
978 {
979 if (bMap1)
980 {
981 ps_rgb_nbox(out, &(mat[i].map[col].rgb), nexty-y);
982 }
983 else
984 {
985 assert(mat2);
986 ps_rgb_nbox(out, &(mat2[i].map[col].rgb), nexty-y);
987 }
988 }
989 else
990 {
991 ps_moverel(out, 0, psr->yboxsize);
992 }
993 y = nexty;
994 bMap1 = bNextMap1;
995 col = nextcol;
996 }
997 }
998 }
999 y0 += box_height(&(mat[i]), psr)+box_dh(psr)+box_dh_top(IS_ONCE, psr);
1000 }
1001
1002 if (psr->X.lineatzero || psr->Y.lineatzero)
1003 {
1004 /* reset y0 for first box */
1005 y0 = dh;
1006 ps_comment(out, "Here starts the zero lines drawing");
1007 draw_zerolines(out, x0, y0, w, nmat, mat, psr);
1008 }
1009
1010 if (elegend != elNone)
1011 {
1012 ps_comment(out, "Now it's legend time!");
1013 ps_linewidth(out, static_cast<int>(psr->linewidth));
1014 if (mat2 == nullptr || elegend != elSecond)
1015 {
1016 bDiscrete = mat[0].bDiscrete;
1017 legend = mat[0].legend;
1018 leg_nmap = nmap1;
1019 leg_map = map1;
1020 }
1021 else
1022 {
1023 bDiscrete = mat2[0].bDiscrete;
1024 legend = mat2[0].legend;
1025 leg_nmap = nmap2;
1026 leg_map = map2;
1027 }
1028 if (bDiscrete)
1029 {
1030 leg_discrete(out, psr->legfontsize, DDD, legend,
1031 psr->legfontsize, psr->legfont, leg_nmap, leg_map);
1032 }
1033 else
1034 {
1035 if (elegend != elBoth)
1036 {
1037 leg_continuous(out, x0+w/2, w/2, DDD, legend,
1038 psr->legfontsize, psr->legfont, leg_nmap, leg_map,
1039 mapoffset);
1040 }
1041 else
1042 {
1043 assert(mat2);
1044 leg_bicontinuous(out, x0+w/2, w, DDD, mat[0].legend, mat2[0].legend,
1045 psr->legfontsize, psr->legfont, nmap1, map1, nmap2, map2);
1046 }
1047 }
1048 ps_comment(out, "Were there, dude");
1049 }
1050
1051 ps_close(out);
1052 }
1053
1054 static void make_axis_labels(int nmat, t_matrix *mat)
1055 {
1056 int i, j;
1057
1058 for (i = 0; (i < nmat); i++)
1059 {
1060 /* Make labels for x axis */
1061 if (mat[i].axis_x == nullptr)
1062 {
1063 snew(mat[i].axis_x, mat[i].nx);
1064 for (j = 0; (j < mat[i].nx); j++)
1065 {
1066 mat[i].axis_x[j] = j;
1067 }
1068 }
1069 /* Make labels for y axis */
1070 if (mat[i].axis_y == nullptr)
1071 {
1072 snew(mat[i].axis_y, mat[i].ny);
1073 for (j = 0; (j < mat[i].ny); j++)
1074 {
1075 mat[i].axis_y[j] = j;
1076 }
1077 }
1078 }
1079 }
1080
1081 static void prune_mat(int nmat, t_matrix *mat, t_matrix *mat2, int skip)
1082 {
1083 int i, x, y, xs, ys;
1084
1085 for (i = 0; i < nmat; i++)
1086 {
1087 fprintf(stderr, "converting %dx%d matrix to %dx%d\n",
1088 mat[i].nx, mat[i].ny,
1089 (mat[i].nx+skip-1)/skip, (mat[i].ny+skip-1)/skip);
1090 /* walk through matrix */
1091 xs = 0;
1092 for (x = 0; (x < mat[i].nx); x++)
1093 {
1094 if (x % skip == 0)
1095 {
1096 mat[i].axis_x[xs] = mat[i].axis_x[x];
1097 if (mat2)
1098 {
1099 mat2[i].axis_x[xs] = mat2[i].axis_x[x];
1100 }
1101 ys = 0;
1102 for (y = 0; (y < mat[i].ny); y++)
1103 {
1104 if (x == 0)
1105 {
1106 mat[i].axis_y[ys] = mat[i].axis_y[y];
1107 if (mat2)
1108 {
1109 mat2[i].axis_y[ys] = mat2[i].axis_y[y];
1110 }
1111 }
1112 if (y % skip == 0)
1113 {
1114 mat[i].matrix[xs][ys] = mat[i].matrix[x][y];
1115 if (mat2)
1116 {
1117 mat2[i].matrix[xs][ys] = mat2[i].matrix[x][y];
1118 }
1119 ys++;
1120 }
1121 }
1122 xs++;
1123 }
1124 }
1125 /* adjust parameters */
1126 mat[i].nx = (mat[i].nx+skip-1)/skip;
1127 mat[i].ny = (mat[i].ny+skip-1)/skip;
1128 if (mat2)
1129 {
1130 mat2[i].nx = (mat2[i].nx+skip-1)/skip;
1131 mat2[i].ny = (mat2[i].ny+skip-1)/skip;
1132 }
1133 }
1134 }
1135
1136 static void zero_lines(int nmat, t_matrix *mat, t_matrix *mat2)
1137 {
1138 int i, x, y, m;
1139 t_matrix *mats;
1140
1141 for (i = 0; i < nmat; i++)
1142 {
1143 for (m = 0; m < (mat2 ? 2 : 1); m++)
1144 {
1145 if (m == 0)
1146 {
1147 mats = mat;
1148 }
1149 else
1150 {
1151 mats = mat2;
1152 }
1153 for (x = 0; x < mats[i].nx-1; x++)
1154 {
1155 if (std::abs(mats[i].axis_x[x+1]) < 1e-5)
1156 {
1157 for (y = 0; y < mats[i].ny; y++)
1158 {
1159 mats[i].matrix[x][y] = 0;
1160 }
1161 }
1162 }
1163 for (y = 0; y < mats[i].ny-1; y++)
1164 {
1165 if (std::abs(mats[i].axis_y[y+1]) < 1e-5)
1166 {
1167 for (x = 0; x < mats[i].nx; x++)
1168 {
1169 mats[i].matrix[x][y] = 0;
1170 }
1171 }
1172 }
1173 }
1174 }
1175 }
1176
1177 static void write_combined_matrix(int ecombine, const char *fn,
1178 int nmat, t_matrix *mat1, t_matrix *mat2,
1179 const real *cmin, const real *cmax)
1180 {
1181 int i, j, k, nlevels;
1182 FILE *out;
1183 real **rmat1, **rmat2;
1184 real rhi, rlo;
1185
1186 out = gmx_ffopen(fn, "w");
1187 for (k = 0; k < nmat; k++)
1188 {
1189 if (mat2[k].nx != mat1[k].nx || mat2[k].ny != mat1[k].ny)
1190 {
1191 gmx_fatal(FARGS, "Size of frame %d in 1st (%dx%d) and 2nd matrix (%dx%d) do"
1192 " not match.\n", k, mat1[k].nx, mat1[k].ny, mat2[k].nx, mat2[k].ny);
1193 }
1194 printf("Combining two %dx%d matrices\n", mat1[k].nx, mat1[k].ny);
1195 rmat1 = matrix2real(&mat1[k], nullptr);
1196 rmat2 = matrix2real(&mat2[k], nullptr);
1197 if (nullptr == rmat1 || nullptr == rmat2)
1198 {
1199 gmx_fatal(FARGS, "Could not extract real data from %s xpm matrices. Note that, e.g.,\n"
1200 "g_rms and g_mdmat provide such data, but not do_dssp.\n",
1201 (nullptr == rmat1 && nullptr == rmat2) ? "both" : "one of the" );
1202 }
1203 rlo = 1e38;
1204 rhi = -1e38;
1205 for (j = 0; j < mat1[k].ny; j++)
1206 {
1207 for (i = 0; i < mat1[k].nx; i++)
1208 {
1209 switch (ecombine)
1210 {
1211 case ecAdd: rmat1[i][j] += rmat2[i][j]; break;
1212 case ecSub: rmat1[i][j] -= rmat2[i][j]; break;
1213 case ecMult: rmat1[i][j] *= rmat2[i][j]; break;
1214 case ecDiv: rmat1[i][j] /= rmat2[i][j]; break;
1215 default:
1216 gmx_fatal(FARGS, "No such combination rule %d for matrices", ecombine);
1217 }
1218 rlo = std::min(rlo, rmat1[i][j]);
1219 rhi = std::max(rhi, rmat1[i][j]);
1220 }
1221 }
1222 if (cmin)
1223 {
1224 rlo = *cmin;
1225 }
1226 if (cmax)
1227 {
1228 rhi = *cmax;
1229 }
1230 nlevels = std::max(mat1[k].nmap, mat2[k].nmap);
1231 if (rhi == rlo)
1232 {
1233 fprintf(stderr,
1234 "combination results in uniform matrix (%g), no output\n", rhi);
1235 }
1236 /*
1237 else if (rlo>=0 || rhi<=0)
1238 write_xpm(out, mat1[k].flags, mat1[k].title, mat1[k].legend,
1239 mat1[k].label_x, mat1[k].label_y,
1240 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1241 rmat1, rlo, rhi, rhi<=0?red:white, rhi<=0?white:blue,
1242 &nlevels);
1243 else
1244 write_xpm3(out, mat2[k].flags, mat1[k].title, mat1[k].legend,
1245 mat1[k].label_x, mat1[k].label_y,
1246 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1247 rmat1, rlo, 0, rhi, red, white, blue, &nlevels);
1248 */
1249 else
1250 {
1251 write_xpm(out, mat1[k].flags, mat1[k].title, mat1[k].legend,
1252 mat1[k].label_x, mat1[k].label_y,
1253 mat1[k].nx, mat1[k].ny, mat1[k].axis_x, mat1[k].axis_y,
1254 rmat1, rlo, rhi, white, black, &nlevels);
1255 }
1256 }
1257 gmx_ffclose(out);
1258 }
1259
1260 static void do_mat(int nmat, t_matrix *mat, t_matrix *mat2,
1261 gmx_bool bFrame, gmx_bool bZeroLine, gmx_bool bDiag, gmx_bool bFirstDiag, gmx_bool bTitle,
1262 gmx_bool bTitleOnce, gmx_bool bYonce, int elegend,
1263 real size, real boxx, real boxy,
1264 const char *epsfile, const char *xpmfile, const char *m2p,
1265 const char *m2pout, int skip, int mapoffset)
1266 {
1267 int i, j, k;
1268
1269 if (mat2)
1270 {
1271 for (k = 0; (k < nmat); k++)
1272 {
1273 if ((mat2[k].nx != mat[k].nx) || (mat2[k].ny != mat[k].ny))
1274 {
1275 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",
1276 k, mat2[k].nx, mat2[k].ny, mat[k].nx, mat[k].ny);
1277 }
1278 for (j = 0; (j < mat[k].ny); j++)
1279 {
1280 for (i = bFirstDiag ? j+1 : j; (i < mat[k].nx); i++)
1281 {
1282 mat[k].matrix[i][j] = mat2[k].matrix[i][j];
1283 }
1284 }
1285 }
1286 }
1287 for (i = 0; (i < nmat); i++)
1288 {
1289 fprintf(stderr, "Matrix %d is %d x %d\n", i, mat[i].nx, mat[i].ny);
1290 }
1291
1292 make_axis_labels(nmat, mat);
1293
1294 if (skip > 1)
1295 {
1296 prune_mat(nmat, mat, mat2, skip);
1297 }
1298
1299 if (bZeroLine)
1300 {
1301 zero_lines(nmat, mat, mat);
1302 }
1303
1304 if (epsfile != nullptr)
1305 {
1306 ps_mat(epsfile, nmat, mat, mat2, bFrame, bDiag, bFirstDiag,
1307 bTitle, bTitleOnce, bYonce, elegend,
1308 size, boxx, boxy, m2p, m2pout, mapoffset);
1309 }
1310 if (xpmfile != nullptr)
1311 {
1312 xpm_mat(xpmfile, nmat, mat, mat2, bDiag, bFirstDiag);
1313 }
1314 }
1315
1316 static void gradient_map(const rvec grad, int nmap, t_mapping map[])
1317 {
1318 int i;
1319 real c;
1320
1321 for (i = 0; i < nmap; i++)
1322 {
1323 c = i/(nmap-1.0);
1324 map[i].rgb.r = 1-c*(1-grad[XX]);
1325 map[i].rgb.g = 1-c*(1-grad[YY]);
1326 map[i].rgb.b = 1-c*(1-grad[ZZ]);
1327 }
1328 }
1329
1330 static void gradient_mat(rvec grad, int nmat, t_matrix mat[])
1331 {
1332 int m;
1333
1334 for (m = 0; m < nmat; m++)
1335 {
1336 gradient_map(grad, mat[m].nmap, mat[m].map);
1337 }
1338 }
1339
1340 static void rainbow_map(gmx_bool bBlue, int nmap, t_mapping map[])
1341 {
1342 int i;
1343 real c, r, g, b;
1344
1345 for (i = 0; i < nmap; i++)
1346 {
1347 c = (map[i].rgb.r + map[i].rgb.g + map[i].rgb.b)/3;
1348 if (c > 1)
1349 {
1350 c = 1;
1351 }
1352 if (bBlue)
1353 {
1354 c = 1 - c;
1355 }
1356 if (c <= 0.25) /* 0-0.25 */
1357 {
1358 r = 0;
1359 g = std::pow(4.0*c, 2.0/3.0);
1360 b = 1;
1361 }
1362 else if (c <= 0.5) /* 0.25-0.5 */
1363 {
1364 r = 0;
1365 g = 1;
1366 b = std::pow(2.0-4.0*c, 2.0/3.0);
1367 }
1368 else if (c <= 0.75) /* 0.5-0.75 */
1369 {
1370 r = std::pow(4.0*c-2.0, 2.0/3.0);
1371 g = 1;
1372 b = 0;
1373 }
1374 else /* 0.75-1 */
1375 {
1376 r = 1;
1377 g = std::pow(4.0-4.0*c, 2.0/3.0);
1378 b = 0;
1379 }
1380 map[i].rgb.r = r;
1381 map[i].rgb.g = g;
1382 map[i].rgb.b = b;
1383 }
1384 }
1385
1386 static void rainbow_mat(gmx_bool bBlue, int nmat, t_matrix mat[])
1387 {
1388 int m;
1389
1390 for (m = 0; m < nmat; m++)
1391 {
1392 rainbow_map(bBlue, mat[m].nmap, mat[m].map);
1393 }
1394 }
1395
1396 int gmx_xpm2ps(int argc, char *argv[])
1397 {
1398 const char *desc[] = {
1399 "[THISMODULE] makes a beautiful color plot of an XPixelMap file.",
1400 "Labels and axis can be displayed, when they are supplied",
1401 "in the correct matrix format.",
1402 "Matrix data may be generated by programs such as [gmx-do_dssp], [gmx-rms] or",
1403 "[gmx-mdmat].[PAR]",
1404 "Parameters are set in the [TT].m2p[tt] file optionally supplied with",
1405 "[TT]-di[tt]. Reasonable defaults are provided. Settings for the [IT]y[it]-axis",
1406 "default to those for the [IT]x[it]-axis. Font names have a defaulting hierarchy:",
1407 "titlefont -> legendfont; titlefont -> (xfont -> yfont -> ytickfont)",
1408 "-> xtickfont, e.g. setting titlefont sets all fonts, setting xfont",
1409 "sets yfont, ytickfont and xtickfont.[PAR]",
1410 "When no [TT].m2p[tt] file is supplied, many settings are taken from",
1411 "command line options. The most important option is [TT]-size[tt],",
1412 "which sets the size of the whole matrix in postscript units.",
1413 "This option can be overridden with the [TT]-bx[tt] and [TT]-by[tt]",
1414 "options (and the corresponding parameters in the [TT].m2p[tt] file),",
1415 "which set the size of a single matrix element.[PAR]",
1416 "With [TT]-f2[tt] a second matrix file can be supplied. Both matrix",
1417 "files will be read simultaneously and the upper left half of the",
1418 "first one ([TT]-f[tt]) is plotted together with the lower right",
1419 "half of the second one ([TT]-f2[tt]). The diagonal will contain",
1420 "values from the matrix file selected with [TT]-diag[tt].",
1421 "Plotting of the diagonal values can be suppressed altogether by",
1422 "setting [TT]-diag[tt] to [TT]none[tt].",
1423 "In this case, a new color map will be generated with",
1424 "a red gradient for negative numbers and a blue for positive.",
1425 "If the color coding and legend labels of both matrices are identical,",
1426 "only one legend will be displayed, else two separate legends are",
1427 "displayed.",
1428 "With [TT]-combine[tt], an alternative operation can be selected",
1429 "to combine the matrices. The output range is automatically set",
1430 "to the actual range of the combined matrix. This can be overridden",
1431 "with [TT]-cmin[tt] and [TT]-cmax[tt].[PAR]",
1432 "[TT]-title[tt] can be set to [TT]none[tt] to suppress the title, or to",
1433 "[TT]ylabel[tt] to show the title in the Y-label position (alongside",
1434 "the [IT]y[it]-axis).[PAR]",
1435 "With the [TT]-rainbow[tt] option, dull grayscale matrices can be turned",
1436 "into attractive color pictures.[PAR]",
1437 "Merged or rainbowed matrices can be written to an XPixelMap file with",
1438 "the [TT]-xpm[tt] option."
1439 };
1440
1441 gmx_output_env_t *oenv;
1442 const char *fn, *epsfile = nullptr, *xpmfile = nullptr;
1443 int i, nmat, nmat2, etitle, elegend, ediag, erainbow, ecombine;
1444 t_matrix *mat = nullptr, *mat2 = nullptr;
1445 gmx_bool bTitle, bTitleOnce, bDiag, bFirstDiag, bGrad;
1446 static gmx_bool bFrame = TRUE, bZeroLine = FALSE, bYonce = FALSE;
1447 static real size = 400, boxx = 0, boxy = 0, cmin = 0, cmax = 0;
1448 static rvec grad = {0, 0, 0};
1449 enum {
1450 etSel, etTop, etOnce, etYlabel, etNone, etNR
1451 };
1452 const char *title[] = { nullptr, "top", "once", "ylabel", "none", nullptr };
1453 /* MUST correspond to enum elXxx as defined at top of file */
1454 const char *legend[] = { nullptr, "both", "first", "second", "none", nullptr };
1455 enum {
1456 edSel, edFirst, edSecond, edNone, edNR
1457 };
1458 const char *diag[] = { nullptr, "first", "second", "none", nullptr };
1459 enum {
1460 erSel, erNo, erBlue, erRed, erNR
1461 };
1462 const char *rainbow[] = { nullptr, "no", "blue", "red", nullptr };
1463 /* MUST correspond to enum ecXxx as defined at top of file */
1464 const char *combine[] = {
1465 nullptr, "halves", "add", "sub", "mult", "div", nullptr
1466 };
1467 static int skip = 1, mapoffset = 0;
1468 t_pargs pa[] = {
1469 { "-frame", FALSE, etBOOL, {&bFrame},
1470 "Display frame, ticks, labels, title and legend" },
1471 { "-title", FALSE, etENUM, {title}, "Show title at" },
1472 { "-yonce", FALSE, etBOOL, {&bYonce}, "Show y-label only once" },
1473 { "-legend", FALSE, etENUM, {legend}, "Show legend" },
1474 { "-diag", FALSE, etENUM, {diag}, "Diagonal" },
1475 { "-size", FALSE, etREAL, {&size},
1476 "Horizontal size of the matrix in ps units" },
1477 { "-bx", FALSE, etREAL, {&boxx},
1478 "Element x-size, overrides [TT]-size[tt] (also y-size when [TT]-by[tt] is not set)" },
1479 { "-by", FALSE, etREAL, {&boxy}, "Element y-size" },
1480 { "-rainbow", FALSE, etENUM, {rainbow},
1481 "Rainbow colors, convert white to" },
1482 { "-gradient", FALSE, etRVEC, {grad},
1483 "Re-scale colormap to a smooth gradient from white {1,1,1} to {r,g,b}" },
1484 { "-skip", FALSE, etINT, {&skip},
1485 "only write out every nr-th row and column" },
1486 { "-zeroline", FALSE, etBOOL, {&bZeroLine},
1487 "insert line in [REF].xpm[ref] matrix where axis label is zero"},
1488 { "-legoffset", FALSE, etINT, {&mapoffset},
1489 "Skip first N colors from [REF].xpm[ref] file for the legend" },
1490 { "-combine", FALSE, etENUM, {combine}, "Combine two matrices" },
1491 { "-cmin", FALSE, etREAL, {&cmin}, "Minimum for combination output" },
1492 { "-cmax", FALSE, etREAL, {&cmax}, "Maximum for combination output" }
1493 };
1494 #define NPA asize(pa)
1495 t_filenm fnm[] = {
1496 { efXPM, "-f", nullptr, ffREAD },
1497 { efXPM, "-f2", "root2", ffOPTRD },
1498 { efM2P, "-di", nullptr, ffLIBOPTRD },
1499 { efM2P, "-do", "out", ffOPTWR },
1500 { efEPS, "-o", nullptr, ffOPTWR },
1501 { efXPM, "-xpm", nullptr, ffOPTWR }
1502 };
1503 #define NFILE asize(fnm)
1504
1505 if (!parse_common_args(&argc, argv, PCA_CAN_VIEW,
1506 NFILE, fnm, NPA, pa,
1507 asize(desc), desc, 0, nullptr, &oenv))
1508 {
1509 return 0;
1510 }
1511
1512 etitle = nenum(title);
1513 elegend = nenum(legend);
1514 ediag = nenum(diag);
1515 erainbow = nenum(rainbow);
1516 ecombine = nenum(combine);
1517 bGrad = opt2parg_bSet("-gradient", NPA, pa);
1518 for (i = 0; i < DIM; i++)
1519 {
1520 if (grad[i] < 0 || grad[i] > 1)
1521 {
1522 gmx_fatal(FARGS, "RGB value %g out of range (0.0-1.0)", grad[i]);
1523 }
1524 }
1525 if (!bFrame)
1526 {
1527 etitle = etNone;
1528 elegend = elNone;
1529 }
1530
1531 epsfile = ftp2fn_null(efEPS, NFILE, fnm);
1532 xpmfile = opt2fn_null("-xpm", NFILE, fnm);
1533 if (epsfile == nullptr && xpmfile == nullptr)
1534 {
1535 if (ecombine != ecHalves)
1536 {
1537 xpmfile = opt2fn("-xpm", NFILE, fnm);
1538 }
1539 else
1540 {
1541 epsfile = ftp2fn(efEPS, NFILE, fnm);
1542 }
1543 }
1544 if (ecombine != ecHalves && epsfile)
1545 {
1546 fprintf(stderr,
1547 "WARNING: can only write result of arithmetic combination "
1548 "of two matrices to .xpm file\n"
1549 " file %s will not be written\n", epsfile);
1550 epsfile = nullptr;
1551 }
1552
1553 bDiag = ediag != edNone;
1554 bFirstDiag = ediag != edSecond;
1555
1556 fn = opt2fn("-f", NFILE, fnm);
1557 nmat = read_xpm_matrix(fn, &mat);
1558 fprintf(stderr, "There %s %d matri%s in %s\n", (nmat > 1) ? "are" : "is", nmat, (nmat > 1) ? "ces" : "x", fn);
1559 fn = opt2fn_null("-f2", NFILE, fnm);
1560 if (fn)
1561 {
1562 nmat2 = read_xpm_matrix(fn, &mat2);
1563 fprintf(stderr, "There %s %d matri%s in %s\n", (nmat2 > 1) ? "are" : "is", nmat2, (nmat2 > 1) ? "ces" : "x", fn);
1564 if (nmat != nmat2)
1565 {
1566 fprintf(stderr, "Different number of matrices, using the smallest number.\n");
1567 nmat = nmat2 = std::min(nmat, nmat2);
1568 }
1569 }
1570 else
1571 {
1572 if (ecombine != ecHalves)
1573 {
1574 fprintf(stderr,
1575 "WARNING: arithmetic matrix combination selected (-combine), "
1576 "but no second matrix (-f2) supplied\n"
1577 " no matrix combination will be performed\n");
1578 }
1579 ecombine = 0;
1580 nmat2 = 0;
1581 }
1582 bTitle = etitle == etTop;
1583 bTitleOnce = etitle == etOnce;
1584 if (etitle == etYlabel)
1585 {
1586 for (i = 0; (i < nmat); i++)
1587 {
1588 std::strcpy(mat[i].label_y, mat[i].title);
1589 if (mat2)
1590 {
1591 std::strcpy(mat2[i].label_y, mat2[i].title);
1592 }
1593 }
1594 }
1595 if (bGrad)
1596 {
1597 gradient_mat(grad, nmat, mat);
1598 if (mat2)
1599 {
1600 gradient_mat(grad, nmat2, mat2);
1601 }
1602 }
1603 if (erainbow != erNo)
1604 {
1605 rainbow_mat(erainbow == erBlue, nmat, mat);
1606 if (mat2)
1607 {
1608 rainbow_mat(erainbow == erBlue, nmat2, mat2);
1609 }
1610 }
1611
1612 if ((mat2 == nullptr) && (elegend != elNone))
1613 {
1614 elegend = elFirst;
1615 }
1616
1617 if (ecombine && ecombine != ecHalves)
1618 {
1619 write_combined_matrix(ecombine, xpmfile, nmat, mat, mat2,
1620 opt2parg_bSet("-cmin", NPA, pa) ? &cmin : nullptr,
1621 opt2parg_bSet("-cmax", NPA, pa) ? &cmax : nullptr);
1622 }
1623 else
1624 {
1625 do_mat(nmat, mat, mat2, bFrame, bZeroLine, bDiag, bFirstDiag,
1626 bTitle, bTitleOnce, bYonce,
1627 elegend, size, boxx, boxy, epsfile, xpmfile,
1628 opt2fn_null("-di", NFILE, fnm), opt2fn_null("-do", NFILE, fnm), skip,
1629 mapoffset);
1630 }
1631
1632 view_all(oenv, NFILE, fnm);
1633
1634 return 0;
1635 }
The ultimate molecular dynamics simulation package