Annotation of OpenXM_contrib/gnuplot/plot3d.c, Revision 1.1
1.1 ! maekawa 1: #ifndef lint
! 2: static char *RCSid = "$Id: plot3d.c,v 1.36 1998/06/18 14:55:14 ddenholm Exp $";
! 3: #endif
! 4:
! 5: /* GNUPLOT - plot3d.c */
! 6:
! 7: /*[
! 8: * Copyright 1986 - 1993, 1998 Thomas Williams, Colin Kelley
! 9: *
! 10: * Permission to use, copy, and distribute this software and its
! 11: * documentation for any purpose with or without fee is hereby granted,
! 12: * provided that the above copyright notice appear in all copies and
! 13: * that both that copyright notice and this permission notice appear
! 14: * in supporting documentation.
! 15: *
! 16: * Permission to modify the software is granted, but not the right to
! 17: * distribute the complete modified source code. Modifications are to
! 18: * be distributed as patches to the released version. Permission to
! 19: * distribute binaries produced by compiling modified sources is granted,
! 20: * provided you
! 21: * 1. distribute the corresponding source modifications from the
! 22: * released version in the form of a patch file along with the binaries,
! 23: * 2. add special version identification to distinguish your version
! 24: * in addition to the base release version number,
! 25: * 3. provide your name and address as the primary contact for the
! 26: * support of your modified version, and
! 27: * 4. retain our contact information in regard to use of the base
! 28: * software.
! 29: * Permission to distribute the released version of the source code along
! 30: * with corresponding source modifications in the form of a patch file is
! 31: * granted with same provisions 2 through 4 for binary distributions.
! 32: *
! 33: * This software is provided "as is" without express or implied warranty
! 34: * to the extent permitted by applicable law.
! 35: ]*/
! 36:
! 37: #include "plot.h"
! 38: #include "setshow.h"
! 39: #include "binary.h"
! 40:
! 41: #ifndef _Windows
! 42: # include "help.h"
! 43: #endif
! 44:
! 45: #ifndef STDOUT
! 46: #define STDOUT 1
! 47: #endif
! 48:
! 49: /* static prototypes */
! 50:
! 51: static void get_3ddata __PROTO((struct surface_points * this_plot));
! 52: static void print_3dtable __PROTO((int pcount));
! 53: static void eval_3dplots __PROTO((void));
! 54: static void grid_nongrid_data __PROTO((struct surface_points * this_plot));
! 55: static void parametric_3dfixup __PROTO((struct surface_points * start_plot, int *plot_num));
! 56:
! 57: /* the curves/surfaces of the plot */
! 58: struct surface_points *first_3dplot = NULL;
! 59: static struct udft_entry plot_func;
! 60:
! 61:
! 62: extern struct udft_entry *dummy_func;
! 63: extern int datatype[];
! 64: extern char timefmt[];
! 65: extern TBOOLEAN is_3d_plot;
! 66: extern int plot_token;
! 67:
! 68: /* in order to support multiple axes, and to
! 69: * simplify ranging in parametric plots, we use
! 70: * arrays to store some things. For 2d plots,
! 71: * elements are y1 = 0 x1 = 1 y2 = 2 x2 = 3
! 72: * for 3d, z = 0, x = 1, y = 2
! 73: * these are given symbolic names in plot.h
! 74: */
! 75:
! 76: extern double min_array[AXIS_ARRAY_SIZE], max_array[AXIS_ARRAY_SIZE];
! 77: extern int auto_array[AXIS_ARRAY_SIZE];
! 78: extern TBOOLEAN log_array[AXIS_ARRAY_SIZE];
! 79: extern double base_array[AXIS_ARRAY_SIZE];
! 80: extern double log_base_array[AXIS_ARRAY_SIZE];
! 81:
! 82: /* some file-wide variables to store which axis we are using */
! 83: static int x_axis, y_axis, z_axis;
! 84:
! 85:
! 86: /* Deleted from setshow.h and renamed */
! 87: extern FILE *gpoutfile;
! 88:
! 89: /* info from datafile module */
! 90: extern int df_datum;
! 91: extern int df_line_number;
! 92: extern int df_no_use_specs;
! 93: extern int df_eof;
! 94: extern int df_timecol[];
! 95: extern TBOOLEAN df_matrix;
! 96:
! 97: #define Inc_c_token if (++c_token >= num_tokens) \
! 98: int_error ("Syntax error", c_token);
! 99:
! 100: /* From plot2d.c */
! 101: extern int reverse_range[];
! 102:
! 103: /*
! 104: * IMHO, code is getting too cluttered with repeated chunks of
! 105: * code. Some macros to simplify, I hope.
! 106: *
! 107: * do { } while(0) is comp.lang.c recommendation for complex macros
! 108: * also means that break can be specified as an action, and it will
! 109: *
! 110: */
! 111:
! 112: /* copy scalar data to arrays
! 113: * optimiser should optimise infinite away
! 114: * dont know we have to support ranges [10:-10] - lets reverse
! 115: * it for now, then fix it at the end.
! 116: */
! 117: #define INIT_ARRAYS(axis, min, max, auto, is_log, base, log_base, infinite) \
! 118: do{if ((auto_array[axis] = auto) == 0 && max<min) {\
! 119: min_array[axis] = max;\
! 120: max_array[axis] = min; /* we will fix later */ \
! 121: } else { \
! 122: min_array[axis] = (infinite && (auto&1)) ? VERYLARGE : min; \
! 123: max_array[axis] = (infinite && (auto&2)) ? -VERYLARGE : max; \
! 124: } \
! 125: log_array[axis] = is_log; base_array[axis] = base; log_base_array[axis] = log_base;\
! 126: }while(0)
! 127:
! 128: /* handle reversed ranges */
! 129: #define CHECK_REVERSE(axis) \
! 130: do{\
! 131: if (auto_array[axis] == 0 && max_array[axis] < min_array[axis]) {\
! 132: double temp = min_array[axis]; min_array[axis] = max_array[axis]; max_array[axis] = temp;\
! 133: reverse_range[axis] = 1; \
! 134: } else reverse_range[axis] = (range_flags[axis]&RANGE_REVERSE); \
! 135: }while(0)
! 136:
! 137:
! 138: /* get optional [min:max] */
! 139: #define LOAD_RANGE(axis) \
! 140: do {\
! 141: if (equals(c_token, "[")) { \
! 142: c_token++; \
! 143: auto_array[axis] = load_range(axis,&min_array[axis], &max_array[axis], auto_array[axis]);\
! 144: if (!equals(c_token, "]"))\
! 145: int_error("']' expected", c_token);\
! 146: c_token++;\
! 147: }\
! 148: } while (0)
! 149:
! 150:
! 151: /* store VALUE or log(VALUE) in STORE, set TYPE as appropriate
! 152: * Do OUT_ACTION or UNDEF_ACTION as appropriate
! 153: * adjust range provided type is INRANGE (ie dont adjust y if x is outrange
! 154: * VALUE must not be same as STORE
! 155: */
! 156:
! 157: #define STORE_WITH_LOG_AND_FIXUP_RANGE(STORE, VALUE, TYPE, AXIS, OUT_ACTION, UNDEF_ACTION)\
! 158: do { if (log_array[AXIS]) { if (VALUE<0.0) {TYPE = UNDEFINED; UNDEF_ACTION; break;} \
! 159: else if (VALUE == 0.0){STORE = -VERYLARGE; TYPE = OUTRANGE; OUT_ACTION; break;} \
! 160: else { STORE = log(VALUE)/log_base_array[AXIS]; } \
! 161: } else STORE = VALUE; \
! 162: if (TYPE != INRANGE) break; /* dont set y range if x is outrange, for example */ \
! 163: if ( VALUE<min_array[AXIS] ) { \
! 164: if (auto_array[AXIS] & 1) min_array[AXIS] = VALUE; else { TYPE = OUTRANGE; OUT_ACTION; break; } \
! 165: } \
! 166: if ( VALUE>max_array[AXIS] ) { \
! 167: if (auto_array[AXIS] & 2) max_array[AXIS] = VALUE; else { TYPE = OUTRANGE; OUT_ACTION; } \
! 168: } \
! 169: } while(0)
! 170:
! 171: /* use this instead empty macro arguments to work around NeXT cpp bug */
! 172: /* if this fails on any system, we might use ((void)0) */
! 173: #define NOOP /* */
! 174:
! 175: /* check range and take logs of min and max if logscale
! 176: * this also restores min and max for ranges like [10:-10]
! 177: */
! 178:
! 179: #ifdef HAVE_STRINGIZE
! 180: # define RANGE_MSG(x) #x " range is less than threshold : see `set zero`"
! 181: #else
! 182: # define RANGE_MSG(x) "x range is less than threshold : see `set zero`"
! 183: #endif
! 184:
! 185: #define FIXUP_RANGE_FOR_LOG(AXIS, WHICH) \
! 186: do { if (reverse_range[AXIS]) { \
! 187: double temp = min_array[AXIS]; \
! 188: min_array[AXIS] = max_array[AXIS]; \
! 189: max_array[AXIS] = temp; \
! 190: }\
! 191: if (log_array[AXIS]) { \
! 192: if (min_array[AXIS] <= 0.0 || max_array[AXIS] <= 0.0) \
! 193: int_error(RANGE_MSG(WHICH), NO_CARET); \
! 194: min_array[AXIS] = log(min_array[AXIS])/log_base_array[AXIS]; \
! 195: max_array[AXIS] = log(max_array[AXIS])/log_base_array[AXIS]; \
! 196: } \
! 197: } while(0)
! 198:
! 199:
! 200:
! 201: /* support for dynamic size of input line */
! 202:
! 203: void plot3drequest()
! 204: /*
! 205: * in the parametric case we would say splot [u= -Pi:Pi] [v= 0:2*Pi] [-1:1]
! 206: * [-1:1] [-1:1] sin(v)*cos(u),sin(v)*cos(u),sin(u) in the non-parametric
! 207: * case we would say only splot [x= -2:2] [y= -5:5] sin(x)*cos(y)
! 208: *
! 209: */
! 210: {
! 211: TBOOLEAN changed;
! 212: int dummy_token0 = -1, dummy_token1 = -1;
! 213:
! 214: is_3d_plot = TRUE;
! 215:
! 216: if (parametric && strcmp(dummy_var[0], "t") == 0) {
! 217: strcpy(dummy_var[0], "u");
! 218: strcpy(dummy_var[1], "v");
! 219: }
! 220: autoscale_lx = autoscale_x;
! 221: autoscale_ly = autoscale_y;
! 222: autoscale_lz = autoscale_z;
! 223:
! 224: if (!term) /* unknown */
! 225: int_error("use 'set term' to set terminal type first", c_token);
! 226:
! 227: if (equals(c_token, "[")) {
! 228: c_token++;
! 229: if (isletter(c_token)) {
! 230: if (equals(c_token + 1, "=")) {
! 231: dummy_token0 = c_token;
! 232: c_token += 2;
! 233: } else {
! 234: /* oops; probably an expression with a variable. */
! 235: /* Parse it as an xmin expression. */
! 236: /* used to be: int_error("'=' expected",c_token); */
! 237: }
! 238: }
! 239: changed = parametric ? load_range(U_AXIS, &umin, &umax, autoscale_lu) : load_range(FIRST_X_AXIS, &xmin, &xmax, autoscale_lx);
! 240: if (!equals(c_token, "]"))
! 241: int_error("']' expected", c_token);
! 242: c_token++;
! 243: /* if (changed) */
! 244: if (parametric)
! 245: /* autoscale_lu = FALSE; */
! 246: autoscale_lu = changed;
! 247: else
! 248: /* autoscale_lx = FALSE; */
! 249: autoscale_lx = changed;
! 250: }
! 251: if (equals(c_token, "[")) {
! 252: c_token++;
! 253: if (isletter(c_token)) {
! 254: if (equals(c_token + 1, "=")) {
! 255: dummy_token1 = c_token;
! 256: c_token += 2;
! 257: } else {
! 258: /* oops; probably an expression with a variable. */
! 259: /* Parse it as an xmin expression. */
! 260: /* used to be: int_error("'=' expected",c_token); */
! 261: }
! 262: }
! 263: changed = parametric ? load_range(V_AXIS, &vmin, &vmax, autoscale_lv) : load_range(FIRST_Y_AXIS, &ymin, &ymax, autoscale_ly);
! 264: if (!equals(c_token, "]"))
! 265: int_error("']' expected", c_token);
! 266: c_token++;
! 267: /* if (changed) */
! 268: if (parametric)
! 269: /* autoscale_lv = FALSE; */
! 270: autoscale_lv = changed;
! 271: else
! 272: /* autoscale_ly = FALSE; */
! 273: autoscale_ly = changed;
! 274: }
! 275: if (parametric) {
! 276: /* set optional x (parametric) or z ranges */
! 277: if (equals(c_token, "[")) {
! 278: c_token++;
! 279: autoscale_lx = load_range(FIRST_X_AXIS, &xmin, &xmax, autoscale_lx);
! 280: if (!equals(c_token, "]"))
! 281: int_error("']' expected", c_token);
! 282: c_token++;
! 283: }
! 284: /* set optional y ranges */
! 285: if (equals(c_token, "[")) {
! 286: c_token++;
! 287: autoscale_ly = load_range(FIRST_Y_AXIS, &ymin, &ymax, autoscale_ly);
! 288: if (!equals(c_token, "]"))
! 289: int_error("']' expected", c_token);
! 290: c_token++;
! 291: }
! 292: } /* parametric */
! 293: if (equals(c_token, "[")) { /* set optional z ranges */
! 294: c_token++;
! 295: autoscale_lz = load_range(FIRST_Z_AXIS, &zmin, &zmax, autoscale_lz);
! 296: if (!equals(c_token, "]"))
! 297: int_error("']' expected", c_token);
! 298: c_token++;
! 299: }
! 300: CHECK_REVERSE(FIRST_X_AXIS);
! 301: CHECK_REVERSE(FIRST_Y_AXIS);
! 302: CHECK_REVERSE(FIRST_Z_AXIS);
! 303:
! 304: /* use the default dummy variable unless changed */
! 305: if (dummy_token0 >= 0)
! 306: copy_str(c_dummy_var[0], dummy_token0, MAX_ID_LEN);
! 307: else
! 308: (void) strcpy(c_dummy_var[0], dummy_var[0]);
! 309:
! 310: if (dummy_token1 >= 0)
! 311: copy_str(c_dummy_var[1], dummy_token1, MAX_ID_LEN);
! 312: else
! 313: (void) strcpy(c_dummy_var[1], dummy_var[1]);
! 314:
! 315: eval_3dplots();
! 316: }
! 317:
! 318:
! 319: static void grid_nongrid_data(this_plot)
! 320: struct surface_points *this_plot;
! 321: {
! 322: int i, j, k;
! 323: double x, y, z, w, dx, dy, xmin, xmax, ymin, ymax;
! 324: struct iso_curve *old_iso_crvs = this_plot->iso_crvs;
! 325: struct iso_curve *icrv, *oicrv, *oicrvs;
! 326:
! 327: /* Compute XY bounding box on the original data. */
! 328: xmin = xmax = old_iso_crvs->points[0].x;
! 329: ymin = ymax = old_iso_crvs->points[0].y;
! 330: for (icrv = old_iso_crvs; icrv != NULL; icrv = icrv->next) {
! 331: struct coordinate GPHUGE *points = icrv->points;
! 332:
! 333: for (i = 0; i < icrv->p_count; i++, points++) {
! 334: if (xmin > points->x)
! 335: xmin = points->x;
! 336: if (xmax < points->x)
! 337: xmax = points->x;
! 338: if (ymin > points->y)
! 339: ymin = points->y;
! 340: if (ymax < points->y)
! 341: ymax = points->y;
! 342: }
! 343: }
! 344:
! 345: dx = (xmax - xmin) / (dgrid3d_col_fineness - 1);
! 346: dy = (ymax - ymin) / (dgrid3d_row_fineness - 1);
! 347:
! 348: /* Create the new grid structure, and compute the low pass filtering from
! 349: * non grid to grid structure.
! 350: */
! 351: this_plot->iso_crvs = NULL;
! 352: this_plot->num_iso_read = dgrid3d_col_fineness;
! 353: this_plot->has_grid_topology = TRUE;
! 354: for (i = 0, x = xmin; i < dgrid3d_col_fineness; i++, x += dx) {
! 355: struct coordinate GPHUGE *points;
! 356:
! 357: icrv = iso_alloc(dgrid3d_row_fineness + 1);
! 358: icrv->p_count = dgrid3d_row_fineness;
! 359: icrv->next = this_plot->iso_crvs;
! 360: this_plot->iso_crvs = icrv;
! 361: points = icrv->points;
! 362:
! 363: for (j = 0, y = ymin; j < dgrid3d_row_fineness; j++, y += dy, points++) {
! 364: z = w = 0.0;
! 365:
! 366: #ifndef BUGGY_DGRID_RANGING /* HBB 981209 */
! 367: /* as soon as ->type is changed to UNDEFINED, break out of
! 368: * two inner loops! */
! 369: points->type = INRANGE;
! 370: #endif
! 371: for (oicrv = old_iso_crvs; oicrv != NULL; oicrv = oicrv->next) {
! 372: struct coordinate GPHUGE *opoints = oicrv->points;
! 373: for (k = 0; k < oicrv->p_count; k++, opoints++) {
! 374: double dist, dist_x = fabs(opoints->x - x), dist_y = fabs(opoints->y - y);
! 375:
! 376: switch (dgrid3d_norm_value) {
! 377: case 1:
! 378: dist = dist_x + dist_y;
! 379: break;
! 380: case 2:
! 381: dist = dist_x * dist_x + dist_y * dist_y;
! 382: break;
! 383: case 4:
! 384: dist = dist_x * dist_x + dist_y * dist_y;
! 385: dist *= dist;
! 386: break;
! 387: case 8:
! 388: dist = dist_x * dist_x + dist_y * dist_y;
! 389: dist *= dist;
! 390: dist *= dist;
! 391: break;
! 392: case 16:
! 393: dist = dist_x * dist_x + dist_y * dist_y;
! 394: dist *= dist;
! 395: dist *= dist;
! 396: dist *= dist;
! 397: break;
! 398: default:
! 399: dist = pow(dist_x, (double) dgrid3d_norm_value) +
! 400: pow(dist_y, (double) dgrid3d_norm_value);
! 401: break;
! 402: }
! 403:
! 404: /* The weight of this point is inverse proportional
! 405: * to the distance.
! 406: */
! 407: if (dist == 0.0) {
! 408: #ifndef BUGGY_DGRID_RANGING
! 409: /* HBB 981209: revised flagging as undefined */
! 410: /* Supporting all those infinities on various
! 411: * platforms becomes tiresome, to say the least :-(
! 412: * Let's just return the first z where this happens,
! 413: * unchanged, and be done with this, period. */
! 414: points->type = UNDEFINED;
! 415: z = opoints->z;
! 416: w = 1.0;
! 417: break; /* out of for (k...) loop */
! 418: #else
! 419: #if !defined(AMIGA_SC_6_1) && !defined(__PUREC__)
! 420: dist = VERYLARGE;
! 421: #else /* !AMIGA_SC_6_1 && !__PUREC__ */
! 422: /* Multiplying VERYLARGE by opoints->z below
! 423: * might yield Inf (i.e. a number that can't
! 424: * be represented on the machine). This will
! 425: * result in points->z being set to NaN. It's
! 426: * better to have a pretty large number that is
! 427: * also on the safe side... The numbers that are
! 428: * read by gnuplot are float values anyway, so
! 429: * they can't be bigger than FLT_MAX. So setting
! 430: * dist to FLT_MAX^2 will make dist pretty large
! 431: * with respect to any value that has been read. */
! 432: dist = ((double) FLT_MAX) * ((double) FLT_MAX);
! 433: #endif /* !AMIGA_SC_6_1 && !__PUREC__ */
! 434: #endif /* BUGGY_DGRID_RANGING */
! 435: } else
! 436: dist = 1.0 / dist;
! 437:
! 438: z += opoints->z * dist;
! 439: w += dist;
! 440: }
! 441: #ifndef BUGGY_DGRID_RANGING
! 442: if (points->type != INRANGE)
! 443: break; /* out of the second-inner loop as well ... */
! 444: #endif
! 445: }
! 446:
! 447: #ifndef BUGGY_DGRID_RANGING
! 448: /* Now that we've escaped the loops safely, we know that we
! 449: * do have a good value in z and w, so we can proceed just as
! 450: * if nothing had happened at all. Nice, isn't it? */
! 451: points->type = INRANGE;
! 452:
! 453: STORE_WITH_LOG_AND_FIXUP_RANGE(points->x, x, points->type, x_axis, NOOP, continue);
! 454: STORE_WITH_LOG_AND_FIXUP_RANGE(points->y, y, points->type, y_axis, NOOP, continue);
! 455: STORE_WITH_LOG_AND_FIXUP_RANGE(points->z, z / w, points->type, z_axis, NOOP, continue);
! 456: #else
! 457: /* HBB 981026: original, short version of this code */
! 458: points->x = x;
! 459: points->y = y;
! 460: points->z = z / w;
! 461: points->type = INRANGE;
! 462: #endif
! 463: }
! 464: }
! 465:
! 466: /* Delete the old non grid data. */
! 467: for (oicrvs = old_iso_crvs; oicrvs != NULL;) {
! 468: oicrv = oicrvs;
! 469: oicrvs = oicrvs->next;
! 470: iso_free(oicrv);
! 471: }
! 472: }
! 473:
! 474: static void get_3ddata(this_plot)
! 475: struct surface_points *this_plot;
! 476: /* this_plot->token is end of datafile spec, before title etc
! 477: * will be moved passed title etc after we return
! 478: */
! 479: {
! 480: int xdatum = 0;
! 481: int ydatum = 0;
! 482: int i, j;
! 483: double v[3];
! 484: int pt_in_iso_crv = 0;
! 485: struct iso_curve *this_iso;
! 486:
! 487: if (mapping3d == MAP3D_CARTESIAN) {
! 488: if (df_no_use_specs == 2)
! 489: int_error("Need 1 or 3 columns for cartesian data", this_plot->token);
! 490: } else {
! 491: if (df_no_use_specs == 1)
! 492: int_error("Need 2 or 3 columns for polar data", this_plot->token);
! 493: }
! 494:
! 495: this_plot->num_iso_read = 0;
! 496: this_plot->has_grid_topology = TRUE;
! 497:
! 498: /* we ought to keep old memory - most likely case
! 499: * is a replot, so it will probably exactly fit into
! 500: * memory already allocated ?
! 501: */
! 502: if (this_plot->iso_crvs != NULL) {
! 503: struct iso_curve *icrv, *icrvs = this_plot->iso_crvs;
! 504:
! 505: while (icrvs) {
! 506: icrv = icrvs;
! 507: icrvs = icrvs->next;
! 508: iso_free(icrv);
! 509: }
! 510: this_plot->iso_crvs = NULL;
! 511: }
! 512: /* data file is already open */
! 513:
! 514: if (df_matrix)
! 515: xdatum = df_3dmatrix(this_plot);
! 516: else {
! 517: /*{{{ read surface from text file */
! 518: struct iso_curve *this_iso = iso_alloc(samples);
! 519: struct coordinate GPHUGE *cp;
! 520: double x, y, z;
! 521:
! 522: while ((j = df_readline(v, 3)) != DF_EOF) {
! 523: if (j == DF_SECOND_BLANK)
! 524: break; /* two blank lines */
! 525: if (j == DF_FIRST_BLANK) {
! 526: /* one blank line */
! 527: if (pt_in_iso_crv == 0) {
! 528: if (xdatum == 0)
! 529: continue;
! 530: pt_in_iso_crv = xdatum;
! 531: }
! 532: if (xdatum > 0) {
! 533: this_iso->p_count = xdatum;
! 534: this_iso->next = this_plot->iso_crvs;
! 535: this_plot->iso_crvs = this_iso;
! 536: this_plot->num_iso_read++;
! 537:
! 538: if (xdatum != pt_in_iso_crv)
! 539: this_plot->has_grid_topology = FALSE;
! 540:
! 541: this_iso = iso_alloc(pt_in_iso_crv);
! 542: xdatum = 0;
! 543: ydatum++;
! 544: }
! 545: continue;
! 546: }
! 547: /* its a data point or undefined */
! 548:
! 549: if (xdatum >= this_iso->p_max) {
! 550: /*
! 551: * overflow about to occur. Extend size of points[] array. We
! 552: * either double the size, or add 1000 points, whichever is a
! 553: * smaller increment. Note i = p_max.
! 554: */
! 555: iso_extend(this_iso,
! 556: xdatum + (xdatum < 1000 ? xdatum : 1000));
! 557: }
! 558: cp = this_iso->points + xdatum;
! 559:
! 560: if (j == DF_UNDEFINED) {
! 561: cp->type = UNDEFINED;
! 562: continue;
! 563: }
! 564: cp->type = INRANGE; /* unless we find out different */
! 565:
! 566: switch (mapping3d) {
! 567: case MAP3D_CARTESIAN:
! 568: switch (j) {
! 569: case 1:
! 570: x = xdatum;
! 571: y = ydatum;
! 572: z = v[0];
! 573: break;
! 574: case 3:
! 575: x = v[0];
! 576: y = v[1];
! 577: z = v[2];
! 578: break;
! 579: default:
! 580: {
! 581: char msg[80];
! 582: sprintf(msg, "Need 1 or 3 columns - line %d", df_line_number);
! 583: int_error(msg, this_plot->token);
! 584: return; /* avoid gcc -Wuninitialised for x,y,z */
! 585: }
! 586: }
! 587: break;
! 588: case MAP3D_SPHERICAL:
! 589: if (j < 2)
! 590: int_error("Need 2 or 3 columns", this_plot->token);
! 591: if (j < 3)
! 592: v[2] = 1; /* default radius */
! 593: if (angles_format == ANGLES_DEGREES) {
! 594: v[0] *= DEG2RAD; /* Convert to radians. */
! 595: v[1] *= DEG2RAD;
! 596: }
! 597: x = v[2] * cos(v[0]) * cos(v[1]);
! 598: y = v[2] * sin(v[0]) * cos(v[1]);
! 599: z = v[2] * sin(v[1]);
! 600: break;
! 601: case MAP3D_CYLINDRICAL:
! 602: if (j < 2)
! 603: int_error("Need 2 or 3 columns", this_plot->token);
! 604: if (j < 3)
! 605: v[2] = 1; /* default radius */
! 606: if (angles_format == ANGLES_DEGREES) {
! 607: v[0] *= DEG2RAD; /* Convert to radians. */
! 608: }
! 609: x = v[2] * cos(v[0]);
! 610: y = v[2] * sin(v[0]);
! 611: z = v[1];
! 612: break;
! 613: default:
! 614: int_error("Internal error : Unknown mapping type", NO_CARET);
! 615: return;
! 616: }
! 617:
! 618: /* adjust for logscales. Set min/max and point types.
! 619: * store in cp
! 620: */
! 621: cp->type = INRANGE;
! 622: /* cannot use continue, as macro is wrapped in a loop.
! 623: * I regard this as correct goto use
! 624: */
! 625: STORE_WITH_LOG_AND_FIXUP_RANGE(cp->x, x, cp->type, x_axis, NOOP, goto come_here_if_undefined);
! 626: STORE_WITH_LOG_AND_FIXUP_RANGE(cp->y, y, cp->type, y_axis, NOOP, goto come_here_if_undefined);
! 627: STORE_WITH_LOG_AND_FIXUP_RANGE(cp->z, z, cp->type, z_axis, NOOP, goto come_here_if_undefined);
! 628:
! 629: /* some may complain, but I regard this as the correct use
! 630: * of goto
! 631: */
! 632: come_here_if_undefined:
! 633: ++xdatum;
! 634: } /* end of whileloop - end of surface */
! 635:
! 636: if (xdatum > 0) {
! 637: this_plot->num_iso_read++; /* Update last iso. */
! 638: this_iso->p_count = xdatum;
! 639:
! 640: this_iso->next = this_plot->iso_crvs;
! 641: this_plot->iso_crvs = this_iso;
! 642:
! 643: if (xdatum != pt_in_iso_crv)
! 644: this_plot->has_grid_topology = FALSE;
! 645:
! 646: } else {
! 647: iso_free(this_iso); /* Free last allocation. */
! 648: }
! 649:
! 650: if (dgrid3d && this_plot->num_iso_read > 0)
! 651: grid_nongrid_data(this_plot);
! 652: /*}}} */
! 653: }
! 654:
! 655: if (this_plot->num_iso_read <= 1)
! 656: this_plot->has_grid_topology = FALSE;
! 657: if (this_plot->has_grid_topology && !hidden3d) {
! 658: struct iso_curve *new_icrvs = NULL;
! 659: int num_new_iso = this_plot->iso_crvs->p_count, len_new_iso = this_plot->num_iso_read;
! 660:
! 661: /* Now we need to set the other direction (pseudo) isolines. */
! 662: for (i = 0; i < num_new_iso; i++) {
! 663: struct iso_curve *new_icrv = iso_alloc(len_new_iso);
! 664:
! 665: new_icrv->p_count = len_new_iso;
! 666:
! 667: for (j = 0, this_iso = this_plot->iso_crvs;
! 668: this_iso != NULL;
! 669: j++, this_iso = this_iso->next) {
! 670: /* copy whole point struct to get type too.
! 671: * wasteful for windows, with padding */
! 672: /* more efficient would be extra pointer to same struct */
! 673: new_icrv->points[j] = this_iso->points[i];
! 674: }
! 675:
! 676: new_icrv->next = new_icrvs;
! 677: new_icrvs = new_icrv;
! 678: }
! 679:
! 680: /* Append the new iso curves after the read ones. */
! 681: for (this_iso = this_plot->iso_crvs;
! 682: this_iso->next != NULL;
! 683: this_iso = this_iso->next);
! 684: this_iso->next = new_icrvs;
! 685: }
! 686: }
! 687:
! 688:
! 689:
! 690: static void print_3dtable(pcount)
! 691: int pcount;
! 692: {
! 693: register struct surface_points *this_plot;
! 694: int i, curve, surface;
! 695: struct iso_curve *icrvs;
! 696: struct coordinate GPHUGE *points;
! 697:
! 698: for (surface = 0, this_plot = first_3dplot; surface < pcount;
! 699: this_plot = this_plot->next_sp, surface++) {
! 700: fprintf(gpoutfile, "\n#Surface %d of %d surfaces\n", surface, pcount);
! 701: icrvs = this_plot->iso_crvs;
! 702: curve = 0;
! 703:
! 704: if (draw_surface) {
! 705: /* only the curves in one direction */
! 706: while (icrvs && curve < this_plot->num_iso_read) {
! 707: fprintf(gpoutfile, "\n#IsoCurve %d, %d points\n#x y z type\n",
! 708: curve, icrvs->p_count);
! 709: for (i = 0, points = icrvs->points; i < icrvs->p_count; i++) {
! 710: fprintf(gpoutfile, "%g %g %g %c\n",
! 711: points[i].x,
! 712: points[i].y,
! 713: points[i].z,
! 714: points[i].type == INRANGE ? 'i'
! 715: : points[i].type == OUTRANGE ? 'o'
! 716: : 'u');
! 717: }
! 718: icrvs = icrvs->next;
! 719: curve++;
! 720: }
! 721: putc('\n', gpoutfile);
! 722: }
! 723: if (draw_contour) {
! 724: int number = 0;
! 725: struct gnuplot_contours *c = this_plot->contours;
! 726: while (c) {
! 727: int count = c->num_pts;
! 728: struct coordinate GPHUGE *p = c->coords;
! 729: if (c->isNewLevel)
! 730: /* dont display count - contour split across chunks */
! 731: /* put # in case user wants to use it for a plot */
! 732: /* double blank line to allow plot ... index ... */
! 733: fprintf(gpoutfile, "\n# Contour %d, label: %s\n", number++, c->label);
! 734: for (; --count >= 0; ++p)
! 735: fprintf(gpoutfile, "%g %g %g\n", p->x, p->y, p->z);
! 736: /* blank line between segments of same contour */
! 737: putc('\n', gpoutfile);
! 738: c = c->next;
! 739: }
! 740: }
! 741: }
! 742: fflush(gpoutfile);
! 743: }
! 744:
! 745:
! 746:
! 747: #define SET_DUMMY_RANGE(AXIS) \
! 748: do{\
! 749: if (parametric || polar) { \
! 750: t_min = tmin; t_max = tmax;\
! 751: } else if (log_array[AXIS]) {\
! 752: if (min_array[AXIS] <= 0.0 || max_array[AXIS] <= 0.0)\
! 753: int_error("x/x2 range must be greater than 0 for log scale!", NO_CARET);\
! 754: t_min = log(min_array[AXIS])/log_base_array[AXIS]; t_max = log(max_array[AXIS])/log_base_array[AXIS];\
! 755: } else {\
! 756: t_min = min_array[AXIS]; t_max = max_array[AXIS];\
! 757: }\
! 758: t_step = (t_max - t_min) / (samples - 1); \
! 759: }while(0)
! 760:
! 761:
! 762: /*
! 763: * This parses the splot command after any range specifications. To support
! 764: * autoscaling on the x/z axis, we want any data files to define the x/y
! 765: * range, then to plot any functions using that range. We thus parse the
! 766: * input twice, once to pick up the data files, and again to pick up the
! 767: * functions. Definitions are processed twice, but that won't hurt.
! 768: * div - okay, it doesn't hurt, but every time an option as added for
! 769: * datafiles, code to parse it has to be added here. Change so that
! 770: * we store starting-token in the plot structure.
! 771: */
! 772: static void eval_3dplots()
! 773: {
! 774: int i, j;
! 775: struct surface_points **tp_3d_ptr;
! 776: int start_token, end_token;
! 777: int begin_token;
! 778: TBOOLEAN some_data_files = FALSE, some_functions = FALSE;
! 779: int plot_num, line_num, point_num, crnt_param = 0; /* 0 = z, 1 = y, 2 = x */
! 780: char *xtitle;
! 781: char *ytitle;
! 782:
! 783: /* Reset first_3dplot. This is usually done at the end of this function.
! 784: * If there is an error within this function, the memory is left allocated,
! 785: * since we cannot call sp_free if the list is incomplete
! 786: */
! 787: first_3dplot = NULL;
! 788:
! 789: /* put stuff into arrays to simplify access */
! 790: INIT_ARRAYS(FIRST_X_AXIS, xmin, xmax, autoscale_lx, is_log_x, base_log_x, log_base_log_x, 0);
! 791: INIT_ARRAYS(FIRST_Y_AXIS, ymin, ymax, autoscale_ly, is_log_y, base_log_y, log_base_log_y, 0);
! 792: INIT_ARRAYS(FIRST_Z_AXIS, zmin, zmax, autoscale_lz, is_log_z, base_log_z, log_base_log_z, 1);
! 793:
! 794: x_axis = FIRST_X_AXIS;
! 795: y_axis = FIRST_Y_AXIS;
! 796: z_axis = FIRST_Z_AXIS;
! 797:
! 798: tp_3d_ptr = &(first_3dplot);
! 799: plot_num = 0;
! 800: line_num = 0; /* default line type */
! 801: point_num = 0; /* default point type */
! 802:
! 803: xtitle = NULL;
! 804: ytitle = NULL;
! 805:
! 806: begin_token = c_token;
! 807:
! 808: /*** First Pass: Read through data files ***/
! 809: /*
! 810: * This pass serves to set the x/yranges and to parse the command, as
! 811: * well as filling in every thing except the function data. That is done
! 812: * after the x/yrange is defined.
! 813: */
! 814: while (TRUE) {
! 815: if (END_OF_COMMAND)
! 816: int_error("function to plt3d expected", c_token);
! 817:
! 818: start_token = c_token;
! 819:
! 820: if (is_definition(c_token)) {
! 821: define();
! 822: } else {
! 823: int specs;
! 824: struct surface_points *this_plot;
! 825:
! 826: if (isstring(c_token)) { /* data file to plot */
! 827:
! 828: /*{{{ data file */
! 829: if (parametric && crnt_param != 0)
! 830: int_error("previous parametric function not fully specified", c_token);
! 831:
! 832: if (!some_data_files) {
! 833: if (autoscale_lx & 1) {
! 834: min_array[FIRST_X_AXIS] = VERYLARGE;
! 835: }
! 836: if (autoscale_lx & 2) {
! 837: max_array[FIRST_X_AXIS] = -VERYLARGE;
! 838: }
! 839: if (autoscale_ly & 1) {
! 840: min_array[FIRST_Y_AXIS] = VERYLARGE;
! 841: }
! 842: if (autoscale_ly & 2) {
! 843: max_array[FIRST_Y_AXIS] = -VERYLARGE;
! 844: }
! 845: some_data_files = TRUE;
! 846: }
! 847: if (*tp_3d_ptr)
! 848: this_plot = *tp_3d_ptr;
! 849: else { /* no memory malloc()'d there yet */
! 850: /* Allocate enough isosamples and samples */
! 851: this_plot = sp_alloc(0, 0, 0, 0);
! 852: *tp_3d_ptr = this_plot;
! 853: }
! 854:
! 855: this_plot->plot_type = DATA3D;
! 856: this_plot->plot_style = data_style;
! 857:
! 858: specs = df_open(3);
! 859: /* parses all datafile-specific modifiers */
! 860: /* we will load the data after parsing title,with,... */
! 861:
! 862: /* for capture to key */
! 863: this_plot->token = end_token = c_token - 1;
! 864:
! 865: /* this_plot->token is temporary, for errors in get_3ddata() */
! 866:
! 867: if (datatype[FIRST_X_AXIS] == TIME) {
! 868: if (specs < 3)
! 869: int_error("Need full using spec for x time data",
! 870: c_token);
! 871: df_timecol[0] = 1;
! 872: }
! 873: if (datatype[FIRST_Y_AXIS] == TIME) {
! 874: if (specs < 3)
! 875: int_error("Need full using spec for y time data",
! 876: c_token);
! 877: df_timecol[1] = 1;
! 878: }
! 879: if (datatype[FIRST_Z_AXIS] == TIME) {
! 880: if (specs < 3)
! 881: df_timecol[0] = 1;
! 882: else
! 883: df_timecol[2] = 1;
! 884: }
! 885: /*}}} */
! 886:
! 887: } else { /* function to plot */
! 888:
! 889: /*{{{ function */
! 890: ++plot_num;
! 891: if (parametric) {
! 892: /* Rotate between x/y/z axes */
! 893: /* +2 same as -1, but beats -ve problem */
! 894: crnt_param = (crnt_param + 2) % 3;
! 895: }
! 896:
! 897: if (*tp_3d_ptr) {
! 898: this_plot = *tp_3d_ptr;
! 899: if (!hidden3d)
! 900: sp_replace(this_plot, samples_1, iso_samples_1,
! 901: samples_2, iso_samples_2);
! 902: else
! 903: sp_replace(this_plot, iso_samples_1, 0,
! 904: 0, iso_samples_2);
! 905:
! 906: } else { /* no memory malloc()'d there yet */
! 907: /* Allocate enough isosamples and samples */
! 908: if (!hidden3d)
! 909: this_plot = sp_alloc(samples_1, iso_samples_1,
! 910: samples_2, iso_samples_2);
! 911: else
! 912: this_plot = sp_alloc(iso_samples_1, 0,
! 913: 0, iso_samples_2);
! 914: *tp_3d_ptr = this_plot;
! 915: }
! 916:
! 917: this_plot->plot_type = FUNC3D;
! 918: this_plot->has_grid_topology = TRUE;
! 919: this_plot->plot_style = func_style;
! 920: this_plot->num_iso_read = iso_samples_2;
! 921: dummy_func = &plot_func;
! 922: plot_func.at = temp_at();
! 923: dummy_func = NULL;
! 924: /* ignore it for now */
! 925: some_functions = TRUE;
! 926: end_token = c_token - 1;
! 927: /*}}} */
! 928:
! 929: } /* end of IS THIS A FILE OR A FUNC block */
! 930:
! 931:
! 932: /*{{{ title */
! 933: if (this_plot->title) {
! 934: free(this_plot->title);
! 935: this_plot->title = NULL;
! 936: }
! 937: if (almost_equals(c_token, "t$itle")) {
! 938: /* if (!isstring(++c_token))
! 939: int_error("Expected title", c_token);
! 940: m_quote_capture(&(this_plot->title), c_token, c_token);
! 941: */
! 942: if (parametric) {
! 943: if (crnt_param != 0)
! 944: int_error("\"title\" allowed only after parametric function fully specified", c_token);
! 945: else {
! 946: if (xtitle != NULL)
! 947: xtitle[0] = NUL; /* Remove default title . */
! 948: if (ytitle != NULL)
! 949: ytitle[0] = NUL; /* Remove default title . */
! 950: }
! 951: }
! 952: if (isstring(++c_token))
! 953: m_quote_capture(&(this_plot->title), c_token, c_token);
! 954: else
! 955: int_error("expecting \"title\" for plot", c_token);
! 956: /* end of new method */
! 957: ++c_token;
! 958: } else if (almost_equals(c_token, "not$itle")) {
! 959: if (xtitle != NULL)
! 960: xtitle[0] = '\0';
! 961: if (ytitle != NULL)
! 962: ytitle[0] = '\0';
! 963: /* this_plot->title = NULL; */
! 964: ++c_token;
! 965: } else {
! 966: m_capture(&(this_plot->title), start_token, end_token);
! 967: if (crnt_param == 2)
! 968: xtitle = this_plot->title;
! 969: else if (crnt_param == 1)
! 970: ytitle = this_plot->title;
! 971: }
! 972: /*}}} */
! 973:
! 974:
! 975: /*{{{ line types, widths, ... */
! 976: this_plot->lp_properties.l_type = line_num;
! 977: this_plot->lp_properties.p_type = point_num;
! 978:
! 979: if (almost_equals(c_token, "w$ith")) {
! 980: this_plot->plot_style = get_style();
! 981: }
! 982: /* pick up line/point specs
! 983: * - point spec allowed if style uses points, ie style&2 != 0
! 984: * - keywords are optional
! 985: */
! 986: LP_PARSE(this_plot->lp_properties, 1, this_plot->plot_style & 2,
! 987: line_num, point_num);
! 988:
! 989: /* allow old-style syntax too - ignore case lt 3 4 for example */
! 990: if (!equals(c_token, ",") && !END_OF_COMMAND) {
! 991: struct value t;
! 992: this_plot->lp_properties.l_type =
! 993: this_plot->lp_properties.p_type = (int) real(const_express(&t)) - 1;
! 994:
! 995: if (!equals(c_token, ",") && !END_OF_COMMAND)
! 996: this_plot->lp_properties.p_type = (int) real(const_express(&t)) - 1;
! 997: }
! 998: if (this_plot->plot_style & 2) /* lines, linesp, ... */
! 999: if (crnt_param == 0)
! 1000: point_num +=
! 1001: 1 + (draw_contour != 0)
! 1002: + (hidden3d != 0);
! 1003:
! 1004: if (crnt_param == 0)
! 1005: line_num += 1 + (draw_contour != 0)
! 1006: + (hidden3d != 0);
! 1007: /*}}} */
! 1008:
! 1009:
! 1010: /* now get the data... having to think hard here...
! 1011: * first time through, we fill in this_plot. For second
! 1012: * surface in file, we have to allocate another surface
! 1013: * struct. BUT we may allocate this store only to
! 1014: * find that it is merely some blank lines at end of file
! 1015: * tp_3d_ptr is still pointing at next field of prev. plot,
! 1016: * before : prev_or_first -> this_plot -> possible_preallocated_store
! 1017: * tp_3d_ptr--^
! 1018: * after : prev_or_first -> first -> second -> last -> possibly_more_store
! 1019: * tp_3d_ptr ----^
! 1020: * if file is empty, tp_3d_ptr is not moved. this_plot continues
! 1021: * to point at allocated storage, but that will be reused later
! 1022: */
! 1023:
! 1024: assert(this_plot == *tp_3d_ptr);
! 1025:
! 1026: if (this_plot->plot_type == DATA3D) {
! 1027: /*{{{ read data */
! 1028: /* remember settings for second surface in file */
! 1029: struct lp_style_type *these_props = &(this_plot->lp_properties);
! 1030: enum PLOT_STYLE this_style = this_plot->plot_style;
! 1031:
! 1032: int this_token = this_plot->token;
! 1033: while (!df_eof) {
! 1034: this_plot = *tp_3d_ptr;
! 1035: assert(this_plot != NULL);
! 1036:
! 1037: /* dont move tp_3d_ptr until we are sure we
! 1038: * have read a surface
! 1039: */
! 1040:
! 1041: /* used by get_3ddata() */
! 1042: this_plot->token = this_token;
! 1043:
! 1044: get_3ddata(this_plot);
! 1045: /* for second pass */
! 1046: this_plot->token = c_token;
! 1047:
! 1048: if (this_plot->num_iso_read == 0)
! 1049: /* probably df_eof, in which case we
! 1050: * will leave loop. if not eof, then
! 1051: * how come we got no surface ? - retry
! 1052: * in neither case do we update tp_3d_ptr
! 1053: */
! 1054: continue;
! 1055:
! 1056: /* okay, we have read a surface */
! 1057: ++plot_num;
! 1058: tp_3d_ptr = &(this_plot->next_sp);
! 1059: if (df_eof)
! 1060: break;
! 1061:
! 1062: /* there might be another surface so allocate
! 1063: * and prepare another surface structure
! 1064: * This does no harm if in fact there are
! 1065: * no more surfaces to read
! 1066: */
! 1067:
! 1068: if ((this_plot = *tp_3d_ptr) != NULL) {
! 1069: if (this_plot->title) {
! 1070: free(this_plot->title);
! 1071: this_plot->title = NULL;
! 1072: }
! 1073: } else {
! 1074: /* Allocate enough isosamples and samples */
! 1075: this_plot = *tp_3d_ptr = sp_alloc(0, 0, 0, 0);
! 1076: }
! 1077:
! 1078: this_plot->plot_type = DATA3D;
! 1079: this_plot->plot_style = this_style;
! 1080: /* Struct copy */
! 1081: this_plot->lp_properties = *these_props;
! 1082: }
! 1083: df_close();
! 1084: /*}}} */
! 1085:
! 1086: } else { /* not a data file */
! 1087: tp_3d_ptr = &(this_plot->next_sp);
! 1088: this_plot->token = c_token; /* store for second pass */
! 1089: }
! 1090:
! 1091: } /* !is_definition() : end of scope of this_plot */
! 1092:
! 1093:
! 1094: if (equals(c_token, ","))
! 1095: c_token++;
! 1096: else
! 1097: break;
! 1098:
! 1099: } /* while(TRUE), ie first pass */
! 1100:
! 1101:
! 1102: if (parametric && crnt_param != 0)
! 1103: int_error("parametric function not fully specified", NO_CARET);
! 1104:
! 1105:
! 1106: /*** Second Pass: Evaluate the functions ***/
! 1107: /*
! 1108: * Everything is defined now, except the function data. We expect no
! 1109: * syntax errors, etc, since the above parsed it all. This makes the code
! 1110: * below simpler. If autoscale_ly, the yrange may still change.
! 1111: * - eh ? - z can still change. x/y/z can change if we are parametric ??
! 1112: */
! 1113:
! 1114: if (some_functions) {
! 1115:
! 1116: /* I've changed the controlled variable in fn plots to u_min etc since
! 1117: * it's easier for me to think parametric - 'normal' plot is after all
! 1118: * a special case. I was confused about x_min being both minimum of
! 1119: * x values found, and starting value for fn plots.
! 1120: */
! 1121: register double u_min, u_max, u_step, v_min, v_max, v_step;
! 1122: double uisodiff, visodiff;
! 1123: struct surface_points *this_plot;
! 1124:
! 1125:
! 1126: if (!parametric) {
! 1127: /*{{{ check ranges */
! 1128: /* give error if xrange badly set from missing datafile error
! 1129: * parametric fn can still set ranges
! 1130: * if there are no fns, we'll report it later as 'nothing to plot'
! 1131: */
! 1132:
! 1133: if (min_array[FIRST_X_AXIS] == VERYLARGE ||
! 1134: max_array[FIRST_X_AXIS] == -VERYLARGE) {
! 1135: int_error("x range is invalid", c_token);
! 1136: }
! 1137: if (min_array[FIRST_Y_AXIS] == VERYLARGE ||
! 1138: max_array[FIRST_Y_AXIS] == -VERYLARGE) {
! 1139: int_error("y range is invalid", c_token);
! 1140: }
! 1141: /* check that xmin -> xmax is not too small */
! 1142: fixup_range(FIRST_X_AXIS, "x");
! 1143: fixup_range(FIRST_Y_AXIS, "y");
! 1144: /*}}} */
! 1145: }
! 1146: if (parametric && !some_data_files) {
! 1147: /*{{{ set ranges */
! 1148: /* parametric fn can still change x/y range */
! 1149: if (autoscale_lx & 1)
! 1150: min_array[FIRST_X_AXIS] = VERYLARGE;
! 1151: if (autoscale_lx & 2)
! 1152: max_array[FIRST_X_AXIS] = -VERYLARGE;
! 1153: if (autoscale_ly & 1)
! 1154: min_array[FIRST_Y_AXIS] = VERYLARGE;
! 1155: if (autoscale_ly & 2)
! 1156: max_array[FIRST_Y_AXIS] = -VERYLARGE;
! 1157: /*}}} */
! 1158: }
! 1159: if (parametric) {
! 1160: u_min = umin;
! 1161: u_max = umax;
! 1162: v_min = vmin;
! 1163: v_max = vmax;
! 1164: } else {
! 1165: /*{{{ figure ranges, taking logs etc into account */
! 1166: if (is_log_x) {
! 1167: if (min_array[FIRST_X_AXIS] <= 0.0 ||
! 1168: max_array[FIRST_X_AXIS] <= 0.0)
! 1169: int_error("x range must be greater than 0 for log scale!",
! 1170: NO_CARET);
! 1171: u_min = log(min_array[FIRST_X_AXIS]) / log_base_log_x;
! 1172: u_max = log(max_array[FIRST_X_AXIS]) / log_base_log_x;
! 1173: } else {
! 1174: u_min = min_array[FIRST_X_AXIS];
! 1175: u_max = max_array[FIRST_X_AXIS];
! 1176: }
! 1177:
! 1178: if (is_log_y) {
! 1179: if (min_array[FIRST_Y_AXIS] <= 0.0 ||
! 1180: max_array[FIRST_Y_AXIS] <= 0.0) {
! 1181: int_error("y range must be greater than 0 for log scale!",
! 1182: NO_CARET);
! 1183: }
! 1184: v_min = log(min_array[FIRST_Y_AXIS]) / log_base_log_y;
! 1185: v_max = log(max_array[FIRST_Y_AXIS]) / log_base_log_y;
! 1186: } else {
! 1187: v_min = min_array[FIRST_Y_AXIS];
! 1188: v_max = max_array[FIRST_Y_AXIS];
! 1189: }
! 1190: /*}}} */
! 1191: }
! 1192:
! 1193:
! 1194: if (samples_1 < 2 || samples_2 < 2 || iso_samples_1 < 2 ||
! 1195: iso_samples_2 < 2) {
! 1196: int_error("samples or iso_samples < 2. Must be at least 2.",
! 1197: NO_CARET);
! 1198: }
! 1199:
! 1200: /* start over */
! 1201: this_plot = first_3dplot;
! 1202: c_token = begin_token;
! 1203:
! 1204: /* why do attributes of this_plot matter ? */
! 1205:
! 1206: if (this_plot && this_plot->has_grid_topology && hidden3d) {
! 1207: u_step = (u_max - u_min) / (iso_samples_1 - 1);
! 1208: v_step = (v_max - v_min) / (iso_samples_2 - 1);
! 1209: } else {
! 1210: u_step = (u_max - u_min) / (samples_1 - 1);
! 1211: v_step = (v_max - v_min) / (samples_2 - 1);
! 1212: }
! 1213:
! 1214: uisodiff = (u_max - u_min) / (iso_samples_1 - 1);
! 1215: visodiff = (v_max - v_min) / (iso_samples_2 - 1);
! 1216:
! 1217:
! 1218: /* Read through functions */
! 1219: while (TRUE) {
! 1220: if (is_definition(c_token)) {
! 1221: define();
! 1222: } else {
! 1223:
! 1224: if (!isstring(c_token)) { /* func to plot */
! 1225: /*{{{ evaluate function */
! 1226: struct iso_curve *this_iso = this_plot->iso_crvs;
! 1227: struct coordinate GPHUGE *points = this_iso->points;
! 1228: int num_sam_to_use, num_iso_to_use;
! 1229:
! 1230: if (parametric)
! 1231: crnt_param = (crnt_param + 2) % 3;
! 1232:
! 1233: dummy_func = &plot_func;
! 1234: plot_func.at = temp_at(); /* reparse function */
! 1235: dummy_func = NULL;
! 1236: num_iso_to_use = iso_samples_2;
! 1237:
! 1238: if (!(this_plot->has_grid_topology && hidden3d))
! 1239: num_sam_to_use = samples_1;
! 1240: else
! 1241: num_sam_to_use = iso_samples_1;
! 1242:
! 1243: for (j = 0; j < num_iso_to_use; j++) {
! 1244: double y = v_min + j * visodiff;
! 1245: /* if (is_log_y) PEM fix logscale y axis */
! 1246: /* y = pow(log_base_log_y,y); 26-Sep-89 */
! 1247: /* parametric => NOT a log quantity (?) */
! 1248: (void) Gcomplex(&plot_func.dummy_values[1],
! 1249: !parametric && is_log_y ? pow(base_log_y, y) : y,
! 1250: 0.0);
! 1251:
! 1252: for (i = 0; i < num_sam_to_use; i++) {
! 1253: double x = u_min + i * u_step;
! 1254: struct value a;
! 1255: double temp;
! 1256:
! 1257: /* if (is_log_x) PEM fix logscale x axis */
! 1258: /* x = pow(base_log_x,x); 26-Sep-89 */
! 1259: /* parametric => NOT a log quantity (?) */
! 1260: (void) Gcomplex(&plot_func.dummy_values[0],
! 1261: !parametric && is_log_x ? pow(base_log_x, x) : x, 0.0);
! 1262:
! 1263: points[i].x = x;
! 1264: points[i].y = y;
! 1265:
! 1266: evaluate_at(plot_func.at, &a);
! 1267:
! 1268: if (undefined || (fabs(imag(&a)) > zero)) {
! 1269: points[i].type = UNDEFINED;
! 1270: continue;
! 1271: }
! 1272: temp = real(&a);
! 1273:
! 1274: points[i].type = INRANGE;
! 1275: STORE_WITH_LOG_AND_FIXUP_RANGE(points[i].z, temp, points[i].type, crnt_param, NOOP, NOOP);
! 1276:
! 1277: }
! 1278: this_iso->p_count = num_sam_to_use;
! 1279: this_iso = this_iso->next;
! 1280: points = this_iso ? this_iso->points : NULL;
! 1281: }
! 1282:
! 1283: if (!(this_plot->has_grid_topology && hidden3d)) {
! 1284: num_iso_to_use = iso_samples_1;
! 1285: num_sam_to_use = samples_2;
! 1286: for (i = 0; i < num_iso_to_use; i++) {
! 1287: double x = u_min + i * uisodiff;
! 1288: /* if (is_log_x) PEM fix logscale x axis */
! 1289: /* x = pow(base_log_x,x); 26-Sep-89 */
! 1290: /* if parametric, no logs involved - 3.6 */
! 1291: (void) Gcomplex(&plot_func.dummy_values[0],
! 1292: (!parametric && is_log_x) ? pow(base_log_x, x) : x, 0.0);
! 1293:
! 1294: for (j = 0; j < num_sam_to_use; j++) {
! 1295: double y = v_min + j * v_step;
! 1296: struct value a;
! 1297: double temp;
! 1298: /* if (is_log_y) PEM fix logscale y axis */
! 1299: /* y = pow(base_log_y,y); 26-Sep-89 */
! 1300: (void) Gcomplex(&plot_func.dummy_values[1],
! 1301: (!parametric && is_log_y) ? pow(base_log_y, y) : y, 0.0);
! 1302:
! 1303: points[j].x = x;
! 1304: points[j].y = y;
! 1305:
! 1306: evaluate_at(plot_func.at, &a);
! 1307:
! 1308: if (undefined || (fabs(imag(&a)) > zero)) {
! 1309: points[j].type = UNDEFINED;
! 1310: continue;
! 1311: }
! 1312: temp = real(&a);
! 1313:
! 1314: points[j].type = INRANGE;
! 1315: STORE_WITH_LOG_AND_FIXUP_RANGE(points[j].z, temp, points[j].type, crnt_param, NOOP, NOOP);
! 1316: }
! 1317: this_iso->p_count = num_sam_to_use;
! 1318: this_iso = this_iso->next;
! 1319: points = this_iso ? this_iso->points : NULL;
! 1320: }
! 1321: }
! 1322: /*}}} */
! 1323: } /* end of ITS A FUNCTION TO PLOT */
! 1324:
! 1325: /* we saved it from first pass */
! 1326: c_token = this_plot->token;
! 1327:
! 1328: /* one data file can make several plots */
! 1329: do
! 1330: this_plot = this_plot->next_sp;
! 1331: while (this_plot && this_plot->token == c_token);
! 1332:
! 1333: } /* !is_definition */
! 1334:
! 1335: if (equals(c_token, ","))
! 1336: c_token++;
! 1337: else
! 1338: break;
! 1339:
! 1340: } /* while(TRUE) */
! 1341:
! 1342:
! 1343: if (parametric) {
! 1344: /* Now actually fix the plot triplets to be single plots. */
! 1345: parametric_3dfixup(first_3dplot, &plot_num);
! 1346: }
! 1347: } /* some functions */
! 1348: /* if first_3dplot is NULL, we have no functions or data at all.
! 1349: * This can happen, if you type "splot x=5", since x=5 is a
! 1350: * variable assignment
! 1351: */
! 1352: if (plot_num == 0 || first_3dplot == NULL) {
! 1353: int_error("no functions or data to plot", c_token);
! 1354: }
! 1355: if (min_array[FIRST_X_AXIS] == VERYLARGE ||
! 1356: max_array[FIRST_X_AXIS] == -VERYLARGE ||
! 1357: min_array[FIRST_Y_AXIS] == VERYLARGE ||
! 1358: max_array[FIRST_Y_AXIS] == -VERYLARGE ||
! 1359: min_array[FIRST_Z_AXIS] == VERYLARGE ||
! 1360: max_array[FIRST_Z_AXIS] == -VERYLARGE)
! 1361: int_error("All points undefined", NO_CARET);
! 1362:
! 1363: fixup_range(FIRST_X_AXIS, "x");
! 1364: fixup_range(FIRST_Y_AXIS, "y");
! 1365: fixup_range(FIRST_Z_AXIS, "z");
! 1366:
! 1367: FIXUP_RANGE_FOR_LOG(FIRST_X_AXIS, x);
! 1368: FIXUP_RANGE_FOR_LOG(FIRST_Y_AXIS, y);
! 1369: FIXUP_RANGE_FOR_LOG(FIRST_Z_AXIS, z);
! 1370:
! 1371: /* last parameter should take plot size into effect...
! 1372: * probably needs to be moved to graph3d.c
! 1373: * in the meantime, a value of 20 gives same behaviour
! 1374: * as 3.5 which will do for the moment
! 1375: */
! 1376:
! 1377: if (xtics)
! 1378: setup_tics(FIRST_X_AXIS, &xticdef, xformat, 20);
! 1379: if (ytics)
! 1380: setup_tics(FIRST_Y_AXIS, &yticdef, yformat, 20);
! 1381: if (ztics)
! 1382: setup_tics(FIRST_Z_AXIS, &zticdef, zformat, 20);
! 1383:
! 1384: #define WRITEBACK(axis,min,max) \
! 1385: if(range_flags[axis]&RANGE_WRITEBACK) \
! 1386: {if (auto_array[axis]&1) min = min_array[axis]; \
! 1387: if (auto_array[axis]&2) max = max_array[axis]; \
! 1388: }
! 1389:
! 1390: WRITEBACK(FIRST_X_AXIS, xmin, xmax);
! 1391: WRITEBACK(FIRST_Y_AXIS, ymin, ymax);
! 1392: WRITEBACK(FIRST_Z_AXIS, zmin, zmax);
! 1393:
! 1394: if (plot_num == 0 || first_3dplot == NULL) {
! 1395: int_error("no functions or data to plot", c_token);
! 1396: }
! 1397: /* Creates contours if contours are to be plotted as well. */
! 1398:
! 1399: if (draw_contour) {
! 1400: struct surface_points *this_plot;
! 1401: for (this_plot = first_3dplot, i = 0;
! 1402: i < plot_num;
! 1403: this_plot = this_plot->next_sp, i++) {
! 1404:
! 1405: if (this_plot->contours) {
! 1406: struct gnuplot_contours *cntrs = this_plot->contours;
! 1407:
! 1408: while (cntrs) {
! 1409: struct gnuplot_contours *cntr = cntrs;
! 1410: cntrs = cntrs->next;
! 1411: free(cntr->coords);
! 1412: free(cntr);
! 1413: }
! 1414: }
! 1415: /* Make sure this one can be contoured. */
! 1416: if (!this_plot->has_grid_topology) {
! 1417: this_plot->contours = NULL;
! 1418: fputs("Notice: cannot contour non grid data!\n", stderr);
! 1419: /* changed from int_error by recommendation of
! 1420: * rkc@xn.ll.mit.edu
! 1421: */
! 1422: } else if (this_plot->plot_type == DATA3D) {
! 1423: this_plot->contours = contour(
! 1424: this_plot->num_iso_read,
! 1425: this_plot->iso_crvs,
! 1426: contour_levels, contour_pts,
! 1427: contour_kind, contour_order,
! 1428: levels_kind, levels_list);
! 1429: } else {
! 1430: this_plot->contours = contour(iso_samples_2,
! 1431: this_plot->iso_crvs,
! 1432: contour_levels, contour_pts,
! 1433: contour_kind, contour_order,
! 1434: levels_kind, levels_list);
! 1435: }
! 1436: }
! 1437: } /* draw_contour */
! 1438: /* perform the plot */
! 1439: if (strcmp(term->name, "table") == 0)
! 1440: print_3dtable(plot_num);
! 1441: else {
! 1442: START_LEAK_CHECK(); /* assert no memory leaks here ! */
! 1443: do_3dplot(first_3dplot, plot_num);
! 1444: END_LEAK_CHECK();
! 1445: }
! 1446:
! 1447: /* if we get here, all went well, so record the line for replot */
! 1448: if (plot_token != -1) {
! 1449: /* note that m_capture also frees the old replot_line */
! 1450: m_capture(&replot_line, plot_token, c_token - 1);
! 1451: plot_token = -1;
! 1452: }
! 1453: sp_free(first_3dplot);
! 1454: first_3dplot = NULL;
! 1455: }
! 1456:
! 1457:
! 1458:
! 1459:
! 1460:
! 1461: static void parametric_3dfixup(start_plot, plot_num)
! 1462: struct surface_points *start_plot;
! 1463: int *plot_num;
! 1464: /*
! 1465: * The hardest part of this routine is collapsing the FUNC plot types in the
! 1466: * list (which are gauranteed to occur in (x,y,z) triplets while preserving
! 1467: * the non-FUNC type plots intact. This means we have to work our way
! 1468: * through various lists. Examples (hand checked):
! 1469: * start_plot:F1->F2->F3->NULL ==> F3->NULL
! 1470: * start_plot:F1->F2->F3->F4->F5->F6->NULL ==> F3->F6->NULL
! 1471: * start_plot:F1->F2->F3->D1->D2->F4->F5->F6->D3->NULL ==>
! 1472: * F3->D1->D2->F6->D3->NULL
! 1473: */
! 1474: {
! 1475: /*
! 1476: * I initialized *free_list with NULL, because my compiler warns some lines
! 1477: * later that it might be uninited. The code however seems to not access that
! 1478: * line in that case, but if I'm right, my change is OK and if not, this is a
! 1479: * serious bug in the code.
! 1480: *
! 1481: * x and y ranges now fixed in eval_3dplots
! 1482: */
! 1483: struct surface_points *xp, *new_list, *free_list = NULL;
! 1484: struct surface_points **last_pointer = &new_list;
! 1485:
! 1486: int i, tlen, surface;
! 1487: char *new_title;
! 1488:
! 1489: /*
! 1490: * Ok, go through all the plots and move FUNC3D types together. Note:
! 1491: * this originally was written to look for a NULL next pointer, but
! 1492: * gnuplot wants to be sticky in grabbing memory and the right number of
! 1493: * items in the plot list is controlled by the plot_num variable.
! 1494: *
! 1495: * Since gnuplot wants to do this sticky business, a free_list of
! 1496: * surface_points is kept and then tagged onto the end of the plot list
! 1497: * as this seems more in the spirit of the original memory behavior than
! 1498: * simply freeing the memory. I'm personally not convinced this sort of
! 1499: * concern is worth it since the time spent computing points seems to
! 1500: * dominate any garbage collecting that might be saved here...
! 1501: */
! 1502: new_list = xp = start_plot;
! 1503: for (surface = 0; surface < *plot_num; surface++) {
! 1504: if (xp->plot_type == FUNC3D) {
! 1505: struct surface_points *yp = xp->next_sp;
! 1506: struct surface_points *zp = yp->next_sp;
! 1507:
! 1508: /* Here's a FUNC3D parametric function defined as three parts.
! 1509: * Go through all the points and assign the x's and y's from xp and
! 1510: * yp to zp. min/max already done
! 1511: */
! 1512: struct iso_curve *xicrvs = xp->iso_crvs;
! 1513: struct iso_curve *yicrvs = yp->iso_crvs;
! 1514: struct iso_curve *zicrvs = zp->iso_crvs;
! 1515:
! 1516: (*plot_num) -= 2;
! 1517:
! 1518: assert(INRANGE < OUTRANGE && OUTRANGE < UNDEFINED);
! 1519:
! 1520: while (zicrvs) {
! 1521: struct coordinate GPHUGE *xpoints = xicrvs->points,
! 1522: GPHUGE * ypoints = yicrvs->points, GPHUGE * zpoints = zicrvs->points;
! 1523: for (i = 0; i < zicrvs->p_count; ++i) {
! 1524: zpoints[i].x = xpoints[i].z;
! 1525: zpoints[i].y = ypoints[i].z;
! 1526: if (zpoints[i].type < xpoints[i].type)
! 1527: zpoints[i].type = xpoints[i].type;
! 1528: if (zpoints[i].type < ypoints[i].type)
! 1529: zpoints[i].type = ypoints[i].type;
! 1530:
! 1531: }
! 1532: xicrvs = xicrvs->next;
! 1533: yicrvs = yicrvs->next;
! 1534: zicrvs = zicrvs->next;
! 1535: }
! 1536:
! 1537: /* Ok, fix up the title to include xp and yp plots. */
! 1538: if (((xp->title && xp->title[0] != '\0') ||
! 1539: (yp->title && yp->title[0] != '\0')) && zp->title) {
! 1540: tlen = (xp->title ? strlen(xp->title) : 0) +
! 1541: (yp->title ? strlen(yp->title) : 0) +
! 1542: (zp->title ? strlen(zp->title) : 0) + 5;
! 1543: new_title = gp_alloc((unsigned long) tlen, "string");
! 1544: new_title[0] = 0;
! 1545: if (xp->title && xp->title[0] != '\0') {
! 1546: strcat(new_title, xp->title);
! 1547: strcat(new_title, ", "); /* + 2 */
! 1548: }
! 1549: if (yp->title && yp->title[0] != '\0') {
! 1550: strcat(new_title, yp->title);
! 1551: strcat(new_title, ", "); /* + 2 */
! 1552: }
! 1553: strcat(new_title, zp->title);
! 1554: free(zp->title);
! 1555: zp->title = new_title;
! 1556: }
! 1557: /* add xp and yp to head of free list */
! 1558: assert(xp->next_sp == yp);
! 1559: yp->next_sp = free_list;
! 1560: free_list = xp;
! 1561:
! 1562: /* add zp to tail of new_list */
! 1563: *last_pointer = zp;
! 1564: last_pointer = &(zp->next_sp);
! 1565: xp = zp->next_sp;
! 1566: } else { /* its a data plot */
! 1567: assert(*last_pointer == xp); /* think this is true ! */
! 1568: last_pointer = &(xp->next_sp);
! 1569: xp = xp->next_sp;
! 1570: }
! 1571: }
! 1572:
! 1573: /* Ok, append free list and write first_plot */
! 1574: *last_pointer = free_list;
! 1575: first_3dplot = new_list;
! 1576: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>