?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Move fileio sources to C++
[gromacs/AngularHB.git] / src / gromacs / fileio / xvgr.cpp
blob1b5a47909aad4411d3e9d90d0a1b19bcb2fc6b80
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, 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 "xvgr.h"
40
41 #include <cctype>
42 #include <cstring>
43
44 #include "gromacs/fileio/gmxfio.h"
45 #include "gromacs/legacyheaders/copyrite.h"
46 #include "gromacs/legacyheaders/oenv.h"
47 #include "gromacs/math/vec.h"
48 #include "gromacs/utility/cstringutil.h"
49 #include "gromacs/utility/exceptions.h"
50 #include "gromacs/utility/fatalerror.h"
51 #include "gromacs/utility/futil.h"
52 #include "gromacs/utility/smalloc.h"
53 #include "gromacs/utility/sysinfo.h"
54
55 gmx_bool output_env_get_print_xvgr_codes(const output_env_t oenv)
56 {
57 int xvg_format;
58
59 xvg_format = output_env_get_xvg_format(oenv);
60
61 return (xvg_format == exvgXMGRACE || xvg_format == exvgXMGR);
62 }
63
64 static char *xvgrstr(const char *gmx, const output_env_t oenv,
65 char *buf, int buflen)
66 {
67 /* Supported greek letter names and corresponding xmgrace/xmgr symbols */
68 const char *sym[] = { "beta", "chi", "delta", "eta", "lambda", "mu", "omega", "phi", "psi", "rho", "theta", NULL };
69 const char symc[] = { 'b', 'c', 'd', 'h', 'l', 'm', 'w', 'f', 'y', 'r', 'q', '\0' };
70 int xvgf;
71 gmx_bool bXVGR;
72 int g, b, i;
73 char c;
74
75 xvgf = output_env_get_xvg_format(oenv);
76 bXVGR = (xvgf == exvgXMGRACE || xvgf == exvgXMGR);
77
78 g = 0;
79 b = 0;
80 while (gmx[g] != '\0')
81 {
82 /* Check with the largest string we have ("lambda"), add one for \0 */
83 if (b + 6 + 1 >= buflen)
84 {
85 gmx_fatal(FARGS, "Output buffer length in xvgstr (%d) too small to process xvg input string '%s'", buflen, gmx);
86 }
87 if (gmx[g] == '\\')
88 {
89 g++;
90 if (gmx[g] == 's')
91 {
92 /* Subscript */
93 if (bXVGR)
94 {
95 buf[b++] = '\\';
96 buf[b++] = 's';
97 }
98 else
99 {
100 buf[b++] = '_';
101 }
102 g++;
103 }
104 else if (gmx[g] == 'S')
105 {
106 /* Superscript */
107 if (bXVGR)
108 {
109 buf[b++] = '\\';
110 buf[b++] = 'S';
111 }
112 else
113 {
114 buf[b++] = '^';
115 }
116 g++;
117 }
118 else if (gmx[g] == 'N')
119 {
120 /* End sub/superscript */
121 if (bXVGR)
122 {
123 buf[b++] = '\\';
124 buf[b++] = 'N';
125 }
126 else
127 {
128 if (gmx[g+1] != ' ')
129 {
130 buf[b++] = ' ';
131 }
132 }
133 g++;
134 }
135 else if (gmx[g] == '4')
136 {
137 /* Backward compatibility for xmgr normal font "\4" */
138 switch (xvgf)
139 {
140 case exvgXMGRACE:
141 sprintf(buf+b, "%s", "\\f{}");
142 break;
143 case exvgXMGR:
144 sprintf(buf+b, "%s", "\\4");
145 break;
146 default:
147 buf[b] = '\0';
148 break;
149 }
150 g++;
151 b = strlen(buf);
152 }
153 else if (gmx[g] == '8')
154 {
155 /* Backward compatibility for xmgr symbol font "\8" */
156 switch (xvgf)
157 {
158 case exvgXMGRACE:
159 sprintf(buf+b, "%s", "\\x");
160 break;
161 case exvgXMGR:
162 sprintf(buf+b, "%s", "\\8");
163 break;
164 default:
165 buf[b] = '\0';
166 break;
167 }
168 g++;
169 b = std::strlen(buf);
170 }
171 else
172 {
173 /* Check for special symbol */
174 i = 0;
175 while (sym[i] != NULL &&
176 gmx_strncasecmp(sym[i], gmx+g, std::strlen(sym[i])) != 0)
177 {
178 i++;
179 }
180 if (sym[i] != NULL)
181 {
182 c = symc[i];
183 if (std::isupper(gmx[g]))
184 {
185 c = std::toupper(c);
186 }
187 switch (xvgf)
188 {
189 case exvgXMGRACE:
190 sprintf(buf+b, "%s%c%s", "\\x", c, "\\f{}");
191 break;
192 case exvgXMGR:
193 sprintf(buf+b, "%s%c%s", "\\8", c, "\\4");
194 break;
195 default:
196 std::strncat(buf+b, gmx+g, std::strlen(sym[i]));
197 b += std::strlen(sym[i]);
198 if (gmx[g+std::strlen(sym[i])] != ' ')
199 {
200 buf[b++] = ' ';
201 }
202 buf[b] = '\0';
203 break;
204 }
205 g += std::strlen(sym[i]);
206 b = std::strlen(buf);
207 }
208 else
209 {
210 /* Unknown escape sequence, this should not happen.
211 * We do not generate a fatal error, since that might
212 * stop otherwise functioning code from working.
213 * Copy the backslash to the output and continue processing.
214 */
215 buf[b++] = '\\';
216 }
217 }
218 }
219 else
220 {
221 buf[b++] = gmx[g++];
222 }
223 }
224
225 buf[b++] = '\0';
226
227 return buf;
228 }
229
230 void xvgr_header(FILE *fp, const char *title, const char *xaxis,
231 const char *yaxis, int exvg_graph_type,
232 const output_env_t oenv)
233 {
234 char pukestr[100], buf[STRLEN];
235
236 if (output_env_get_print_xvgr_codes(oenv))
237 {
238 gmx_format_current_time(buf, STRLEN);
239 fprintf(fp, "# This file was created %s", buf);
240 try
241 {
242 gmx::BinaryInformationSettings settings;
243 settings.generatedByHeader(true);
244 settings.linePrefix("# ");
245 gmx::printBinaryInformation(fp, output_env_get_program_context(oenv),
246 settings);
247 }
248 GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR;
249 fprintf(fp, "# %s is part of G R O M A C S:\n#\n",
250 output_env_get_program_display_name(oenv));
251 bromacs(pukestr, 99);
252 fprintf(fp, "# %s\n#\n", pukestr);
253 fprintf(fp, "@ title \"%s\"\n", xvgrstr(title, oenv, buf, STRLEN));
254 fprintf(fp, "@ xaxis label \"%s\"\n",
255 xvgrstr(xaxis, oenv, buf, STRLEN));
256 fprintf(fp, "@ yaxis label \"%s\"\n",
257 xvgrstr(yaxis, oenv, buf, STRLEN));
258 switch (exvg_graph_type)
259 {
260 case exvggtXNY:
261 if (output_env_get_xvg_format(oenv) == exvgXMGR)
262 {
263 fprintf(fp, "@TYPE nxy\n");
264 }
265 else
266 {
267 fprintf(fp, "@TYPE xy\n");
268 }
269 break;
270 case exvggtXYDY:
271 fprintf(fp, "@TYPE xydy\n");
272 break;
273 case exvggtXYDYDY:
274 fprintf(fp, "@TYPE xydydy\n");
275 break;
276 }
277 }
278 }
279
280 FILE *xvgropen_type(const char *fn, const char *title, const char *xaxis,
281 const char *yaxis, int exvg_graph_type,
282 const output_env_t oenv)
283 {
284 FILE *fp;
285
286 fp = gmx_fio_fopen(fn, "w");
287
288 xvgr_header(fp, title, xaxis, yaxis, exvg_graph_type, oenv);
289
290 return fp;
291 }
292
293 FILE *xvgropen(const char *fn, const char *title, const char *xaxis,
294 const char *yaxis, const output_env_t oenv)
295 {
296 return xvgropen_type(fn, title, xaxis, yaxis, exvggtXNY, oenv);
297 }
298
299 void
300 xvgrclose(FILE *fp)
301 {
302 gmx_fio_fclose(fp);
303 }
304
305 void xvgr_subtitle(FILE *out, const char *subtitle, const output_env_t oenv)
306 {
307 char buf[STRLEN];
308
309 if (output_env_get_print_xvgr_codes(oenv))
310 {
311 fprintf(out, "@ subtitle \"%s\"\n", xvgrstr(subtitle, oenv, buf, STRLEN));
312 }
313 }
314
315 void xvgr_view(FILE *out, real xmin, real ymin, real xmax, real ymax,
316 const output_env_t oenv)
317 {
318 if (output_env_get_print_xvgr_codes(oenv))
319 {
320 fprintf(out, "@ view %g, %g, %g, %g\n", xmin, ymin, xmax, ymax);
321 }
322 }
323
324 void xvgr_world(FILE *out, real xmin, real ymin, real xmax, real ymax,
325 const output_env_t oenv)
326 {
327 if (output_env_get_print_xvgr_codes(oenv))
328 {
329 fprintf(out, "@ world xmin %g\n"
330 "@ world ymin %g\n"
331 "@ world xmax %g\n"
332 "@ world ymax %g\n", xmin, ymin, xmax, ymax);
333 }
334 }
335
336 void xvgr_legend(FILE *out, int nsets, const char** setname,
337 const output_env_t oenv)
338 {
339 int i;
340 char buf[STRLEN];
341
342 if (output_env_get_print_xvgr_codes(oenv))
343 {
344 xvgr_view(out, 0.15, 0.15, 0.75, 0.85, oenv);
345 fprintf(out, "@ legend on\n");
346 fprintf(out, "@ legend box on\n");
347 fprintf(out, "@ legend loctype view\n");
348 fprintf(out, "@ legend %g, %g\n", 0.78, 0.8);
349 fprintf(out, "@ legend length %d\n", 2);
350 for (i = 0; (i < nsets); i++)
351 {
352 if (setname[i])
353 {
354 if (output_env_get_xvg_format(oenv) == exvgXMGR)
355 {
356 fprintf(out, "@ legend string %d \"%s\"\n",
357 i, xvgrstr(setname[i], oenv, buf, STRLEN));
358 }
359 else
360 {
361 fprintf(out, "@ s%d legend \"%s\"\n",
362 i, xvgrstr(setname[i], oenv, buf, STRLEN));
363 }
364 }
365 }
366 }
367 }
368
369 void xvgr_new_dataset(FILE *out, int nr_first, int nsets,
370 const char **setname,
371 const output_env_t oenv)
372 {
373 int i;
374 char buf[STRLEN];
375
376 if (output_env_get_print_xvgr_codes(oenv))
377 {
378 fprintf(out, "@\n");
379 for (i = 0; (i < nsets); i++)
380 {
381 if (setname[i])
382 {
383 if (output_env_get_xvg_format(oenv) == exvgXMGR)
384 {
385 fprintf(out, "@ legend string %d \"%s\"\n",
386 i+nr_first, xvgrstr(setname[i], oenv, buf, STRLEN));
387 }
388 else
389 {
390 fprintf(out, "@ s%d legend \"%s\"\n",
391 i+nr_first, xvgrstr(setname[i], oenv, buf, STRLEN));
392 }
393 }
394 }
395 }
396 else
397 {
398 fprintf(out, "\n");
399 }
400 }
401
402 void xvgr_line_props(FILE *out, int NrSet, int LineStyle, int LineColor,
403 const output_env_t oenv)
404 {
405 if (output_env_get_print_xvgr_codes(oenv))
406 {
407 fprintf(out, "@ with g0\n");
408 fprintf(out, "@ s%d linestyle %d\n", NrSet, LineStyle);
409 fprintf(out, "@ s%d color %d\n", NrSet, LineColor);
410 }
411 }
412
413 static const char *LocTypeStr[] = { "view", "world" };
414 static const char *BoxFillStr[] = { "none", "color", "pattern" };
415
416 void xvgr_box(FILE *out,
417 int LocType,
418 real xmin, real ymin, real xmax, real ymax,
419 int LineStyle, int LineWidth, int LineColor,
420 int BoxFill, int BoxColor, int BoxPattern, const output_env_t oenv)
421 {
422 if (output_env_get_print_xvgr_codes(oenv))
423 {
424 fprintf(out, "@with box\n");
425 fprintf(out, "@ box on\n");
426 fprintf(out, "@ box loctype %s\n", LocTypeStr[LocType]);
427 fprintf(out, "@ box %g, %g, %g, %g\n", xmin, ymin, xmax, ymax);
428 fprintf(out, "@ box linestyle %d\n", LineStyle);
429 fprintf(out, "@ box linewidth %d\n", LineWidth);
430 fprintf(out, "@ box color %d\n", LineColor);
431 fprintf(out, "@ box fill %s\n", BoxFillStr[BoxFill]);
432 fprintf(out, "@ box fill color %d\n", BoxColor);
433 fprintf(out, "@ box fill pattern %d\n", BoxPattern);
434 fprintf(out, "@box def\n");
435 }
436 }
437
438 /* reads a line into ptr, adjusting len and renewing ptr if neccesary */
439 static char *fgets3(FILE *fp, char **ptr, int *len, int maxlen)
440 {
441 int len_remaining = *len; /* remaining amount of allocated bytes in buf */
442 int curp = 0; /* current position in buf to read into */
443
444 do
445 {
446 if (len_remaining < 2)
447 {
448 if (*len + STRLEN < maxlen)
449 {
450 /* This line is longer than len characters, let's increase len! */
451 *len += STRLEN;
452 len_remaining += STRLEN;
453 srenew(*ptr, *len);
454 }
455 else
456 {
457 /*something is wrong, we'll just keep reading and return NULL*/
458 len_remaining = STRLEN;
459 curp = 0;
460 }
461 }
462 if (fgets(*ptr + curp, len_remaining, fp) == NULL)
463 {
464 /* if last line, skip */
465 return NULL;
466 }
467 curp += len_remaining-1; /* overwrite the nul char in next iteration */
468 len_remaining = 1;
469 }
470 while ((std::strchr(*ptr, '\n') == NULL) && (!feof(fp)));
471
472 if (*len + STRLEN >= maxlen)
473 {
474 return NULL; /* this line was too long */
475 }
476
477 if (feof(fp))
478 {
479 /* We reached EOF before '\n', skip this last line. */
480 return NULL;
481 }
482 {
483 /* now remove newline */
484 int slen = std::strlen(*ptr);
485 if ((*ptr)[slen-1] == '\n')
486 {
487 (*ptr)[slen-1] = '\0';
488 }
489 }
490
491 return *ptr;
492 }
493
494 static int wordcount(char *ptr)
495 {
496 int i, n = 0, is[2];
497 int cur = 0;
498 #define prev (1-cur)
499
500 if (NULL != ptr)
501 {
502 for (i = 0; (ptr[i] != '\0'); i++)
503 {
504 is[cur] = std::isspace(ptr[i]);
505 if ((0 == i) && !is[cur])
506 {
507 n++;
508 }
509 else if ((i > 0) && (!is[cur] && is[prev]))
510 {
511 n++;
512 }
513 cur = prev;
514 }
515 }
516 return n;
517 }
518
519 static char *read_xvgr_string(const char *line)
520 {
521 const char *ptr0, *ptr1;
522 char *str;
523
524 ptr0 = std::strchr(line, '"');
525 if (ptr0 != NULL)
526 {
527 ptr0++;
528 ptr1 = std::strchr(ptr0, '"');
529 if (ptr1 != NULL)
530 {
531 str = gmx_strdup(ptr0);
532 str[ptr1-ptr0] = '\0';
533 }
534 else
535 {
536 str = gmx_strdup("");
537 }
538 }
539 else
540 {
541 str = gmx_strdup("");
542 }
543
544 return str;
545 }
546
547 int read_xvg_legend(const char *fn, double ***y, int *ny,
548 char **subtitle, char ***legend)
549 {
550 FILE *fp;
551 char *ptr;
552 char *base = NULL;
553 char *fmt = NULL;
554 int k, line = 0, nny, nx, maxx, rval, legend_nalloc, set, nchar;
555 double lf;
556 double **yy = NULL;
557 char *tmpbuf;
558 int len = STRLEN;
559 *ny = 0;
560 nny = 0;
561 nx = 0;
562 maxx = 0;
563 fp = gmx_fio_fopen(fn, "r");
564
565 snew(tmpbuf, len);
566 if (subtitle != NULL)
567 {
568 *subtitle = NULL;
569 }
570 legend_nalloc = 0;
571 if (legend != NULL)
572 {
573 *legend = NULL;
574 }
575
576 while ((ptr = fgets3(fp, &tmpbuf, &len, 10*STRLEN)) != NULL && ptr[0] != '&')
577 {
578 line++;
579 trim(ptr);
580 if (ptr[0] == '@')
581 {
582 if (legend != NULL)
583 {
584 ptr++;
585 trim(ptr);
586 set = -1;
587 if (std::strncmp(ptr, "subtitle", 8) == 0)
588 {
589 ptr += 8;
590 if (subtitle != NULL)
591 {
592 *subtitle = read_xvgr_string(ptr);
593 }
594 }
595 else if (std::strncmp(ptr, "legend string", 13) == 0)
596 {
597 ptr += 13;
598 sscanf(ptr, "%d%n", &set, &nchar);
599 ptr += nchar;
600 }
601 else if (ptr[0] == 's')
602 {
603 ptr++;
604 sscanf(ptr, "%d%n", &set, &nchar);
605 ptr += nchar;
606 trim(ptr);
607 if (std::strncmp(ptr, "legend", 6) == 0)
608 {
609 ptr += 6;
610 }
611 else
612 {
613 set = -1;
614 }
615 }
616 if (set >= 0)
617 {
618 if (set >= legend_nalloc)
619 {
620 legend_nalloc = set + 1;
621 srenew(*legend, legend_nalloc);
622 (*legend)[set] = read_xvgr_string(ptr);
623 }
624 }
625 }
626 }
627 else if (ptr[0] != '#')
628 {
629 if (nny == 0)
630 {
631 (*ny) = nny = wordcount(ptr);
632 /* fprintf(stderr,"There are %d columns in your file\n",nny);*/
633 if (nny == 0)
634 {
635 return 0;
636 }
637 snew(yy, nny);
638 snew(fmt, 3*nny+1);
639 snew(base, 3*nny+1);
640 }
641 /* Allocate column space */
642 if (nx >= maxx)
643 {
644 maxx += 1024;
645 for (k = 0; (k < nny); k++)
646 {
647 srenew(yy[k], maxx);
648 }
649 }
650 /* Initiate format string */
651 fmt[0] = '\0';
652 base[0] = '\0';
653
654 /* fprintf(stderr,"ptr='%s'\n",ptr);*/
655 for (k = 0; (k < nny); k++)
656 {
657 std::strcpy(fmt, base);
658 std::strcat(fmt, "%lf");
659 rval = sscanf(ptr, fmt, &lf);
660 /* fprintf(stderr,"rval = %d\n",rval);*/
661 if ((rval == EOF) || (rval == 0))
662 {
663 break;
664 }
665 yy[k][nx] = lf;
666 srenew(fmt, 3*(nny+1)+1);
667 srenew(base, 3*nny+1);
668 std::strcat(base, "%*s");
669 }
670 if (k != nny)
671 {
672 fprintf(stderr, "Only %d columns on line %d in file %s\n",
673 k, line, fn);
674 for (; (k < nny); k++)
675 {
676 yy[k][nx] = 0.0;
677 }
678 }
679 nx++;
680 }
681 }
682 gmx_fio_fclose(fp);
683
684 *y = yy;
685 sfree(tmpbuf);
686 sfree(base);
687 sfree(fmt);
688
689 if (legend_nalloc > 0)
690 {
691 if (*ny - 1 > legend_nalloc)
692 {
693 srenew(*legend, *ny-1);
694 for (set = legend_nalloc; set < *ny-1; set++)
695 {
696 (*legend)[set] = NULL;
697 }
698 }
699 }
700
701 return nx;
702 }
703
704 int read_xvg(const char *fn, double ***y, int *ny)
705 {
706 return read_xvg_legend(fn, y, ny, NULL, NULL);
707 }
708
709 void write_xvg(const char *fn, const char *title, int nx, int ny, real **y,
710 const char **leg, const output_env_t oenv)
711 {
712 FILE *fp;
713 int i, j;
714
715 fp = xvgropen(fn, title, "X", "Y", oenv);
716 if (leg)
717 {
718 xvgr_legend(fp, ny-1, leg, oenv);
719 }
720 for (i = 0; (i < nx); i++)
721 {
722 for (j = 0; (j < ny); j++)
723 {
724 fprintf(fp, " %12.5e", y[j][i]);
725 }
726 fprintf(fp, "\n");
727 }
728 xvgrclose(fp);
729 }
730
731 real **read_xvg_time(const char *fn,
732 gmx_bool bHaveT, gmx_bool bTB, real tb, gmx_bool bTE, real te,
733 int nsets_in, int *nset, int *nval, real *dt, real **t)
734 {
735 FILE *fp;
736 #define MAXLINELEN 16384
737 char line0[MAXLINELEN];
738 char *line;
739 int t_nalloc, *val_nalloc, a, narg, n, sin, set, nchar;
740 double dbl;
741 gmx_bool bEndOfSet, bTimeInRange, bFirstLine = TRUE;
742 real **val;
743
744 t_nalloc = 0;
745 *t = NULL;
746 val = NULL;
747 val_nalloc = NULL;
748 *nset = 0;
749 *dt = 0;
750 fp = gmx_fio_fopen(fn, "r");
751 for (sin = 0; sin < nsets_in; sin++)
752 {
753 if (nsets_in == 1)
754 {
755 narg = 0;
756 }
757 else
758 {
759 narg = bHaveT ? 2 : 1;
760 }
761 n = 0;
762 bEndOfSet = FALSE;
763 while (!bEndOfSet && fgets(line0, MAXLINELEN, fp))
764 {
765 line = line0;
766 /* Remove whitespace */
767 while (line[0] == ' ' || line[0] == '\t')
768 {
769 line++;
770 }
771 bEndOfSet = (line[0] == '&');
772 if (line[0] != '#' && line[0] != '@' && line[0] != '\n' && !bEndOfSet)
773 {
774 if (bFirstLine && bHaveT)
775 {
776 /* Check the first line that should contain data */
777 a = sscanf(line, "%lf%lf", &dbl, &dbl);
778 if (a == 0)
779 {
780 gmx_fatal(FARGS, "Expected a number in %s on line:\n%s", fn, line0);
781 }
782 else if (a == 1)
783 {
784 fprintf(stderr, "Found only 1 number on line, "
785 "assuming no time is present.\n");
786 bHaveT = FALSE;
787 if (nsets_in > 1)
788 {
789 narg = 1;
790 }
791 }
792 }
793
794 a = 0;
795 bTimeInRange = TRUE;
796 while ((a < narg || (nsets_in == 1 && n == 0)) &&
797 sscanf(line, "%lf%n", &dbl, &nchar) == 1 && bTimeInRange)
798 {
799 /* Use set=-1 as the time "set" */
800 if (sin)
801 {
802 if (!bHaveT || (a > 0))
803 {
804 set = sin;
805 }
806 else
807 {
808 set = -1;
809 }
810 }
811 else
812 {
813 if (!bHaveT)
814 {
815 set = a;
816 }
817 else
818 {
819 set = a-1;
820 }
821 }
822 if (set == -1 && ((bTB && dbl < tb) || (bTE && dbl > te)))
823 {
824 bTimeInRange = FALSE;
825 }
826
827 if (bTimeInRange)
828 {
829 if (n == 0)
830 {
831 if (nsets_in == 1)
832 {
833 narg++;
834 }
835 if (set >= 0)
836 {
837 *nset = set+1;
838 srenew(val, *nset);
839 srenew(val_nalloc, *nset);
840 val_nalloc[set] = 0;
841 val[set] = NULL;
842 }
843 }
844 if (set == -1)
845 {
846 if (sin == 0)
847 {
848 if (n >= t_nalloc)
849 {
850 t_nalloc = over_alloc_small(n);
851 srenew(*t, t_nalloc);
852 }
853 (*t)[n] = dbl;
854 }
855 /* else we should check the time of the next sets with set 0 */
856 }
857 else
858 {
859 if (n >= val_nalloc[set])
860 {
861 val_nalloc[set] = over_alloc_small(n);
862 srenew(val[set], val_nalloc[set]);
863 }
864 val[set][n] = (real)dbl;
865 }
866 }
867 a++;
868 line += nchar;
869 }
870 if (line0[strlen(line0)-1] != '\n')
871 {
872 fprintf(stderr, "File %s does not end with a newline, ignoring the last line\n", fn);
873 }
874 else if (bTimeInRange)
875 {
876 if (a == 0)
877 {
878 fprintf(stderr, "Ignoring invalid line in %s:\n%s", fn, line0);
879 }
880 else
881 {
882 if (a != narg)
883 {
884 fprintf(stderr, "Invalid line in %s:\n%s"
885 "Using zeros for the last %d sets\n",
886 fn, line0, narg-a);
887 }
888 n++;
889 }
890 }
891 if (a > 0)
892 {
893 bFirstLine = FALSE;
894 }
895 }
896 }
897 if (sin == 0)
898 {
899 *nval = n;
900 if (!bHaveT)
901 {
902 snew(*t, n);
903 for (a = 0; a < n; a++)
904 {
905 (*t)[a] = a;
906 }
907 }
908 if (n > 1)
909 {
910 *dt = (real)((*t)[n-1]-(*t)[0])/(n-1.0);
911 }
912 else
913 {
914 *dt = 1;
915 }
916 }
917 else
918 {
919 if (n < *nval)
920 {
921 fprintf(stderr, "Set %d is shorter (%d) than the previous set (%d)\n",
922 sin+1, n, *nval);
923 *nval = n;
924 fprintf(stderr, "Will use only the first %d points of every set\n",
925 *nval);
926 }
927 }
928 }
929 gmx_fio_fclose(fp);
930
931 sfree(val_nalloc);
932
933 return val;
934 }