xstarplot.C Source File

00001 /*
00002  *  xstarplot.C:  plot an N-body system in an X environment.
00003  *
00004  *.............................................................................
00005  *
00006  *    version 1:  May 1989   Piet Hut            email: piet@iassns.bitnet
00007  *                           Institute for Advanced Study, Princeton, NJ, USA
00008  *    version 2:  Dec 1992   Piet Hut      adapted to the new C++-based Starlab
00009  *    version 3:  Dec 1992   Steve McMillan      email: steve@zonker.drexel.edu
00010  *                           Drexel University, Philadelphia, PA, USA
00011  *                           Converted original ASCIIplot to X --- FIRST DRAFT!
00012  *    version 4:  Jan 1994   Biao Lu             email: biao@eagle.drexel.edu
00013  *                           Drexel University, Philadelphia, PA, USA
00014  *                           Use new X interface for movie.
00015  *            Aug/Sep 1994   Cleaned up and expanded, SMcM
00016  *                Sep 1995   Extended to sdyn3 and nonzero radii, SMcM
00017  *
00018  *.............................................................................
00019  *
00020  *  non-local functions: 
00021  *    xstarplot
00022  *
00023  *.............................................................................
00024  *
00025  *  uses:
00026  *    liblux.a (X interface)
00027  *    The X graphics library
00028  *
00029  *.............................................................................
00030  */
00031 
00032 // NOTE: On some systems, it appears that overflow can occur, causing the
00033 //       win.lnx and win.lny entries to be overwritten and leading to
00034 //       spurious "log10" errors in draw.c (in lux_set_item, specifically).
00035 //       This isn't fatal, but it should be fixed someday...
00036 
00037 #include "sdyn3.h"
00038 #include "xstarplot.h"
00039 #include "dyn_util.h"
00040 
00041 //-------------------------------------------------------------------------
00042 
00043 // These shouldn't really be global, but keep them this way for now...
00044 
00045 unsigned long  win, dia, instr, colwin;
00046 unsigned long c_energy[10], c_index[N_COLORS+1];
00047 int      init_status = 0;
00048 
00049 // The plotting box is determined completely by origin and lmax3d.
00050 // Its projection onto the viewing plane is xmin, xmax, ymin, ymax
00051 
00052 float  xmin, xmax, ymin, ymax, base_point_size, lmax3d, origin[3];
00053 
00054 int   win_size;
00055 int   point_scale_mode = 0;
00056 int   kx = 1, ky = 2, kproj = 3;
00057 int   origin_star = -1;
00058 int   nodes = 0, links = 0, root = 0;
00059 float theta = 0.33, costheta = cos(0.33), sintheta = sin(0.33), dtheta = 0.03;
00060 float phi = 0.33, cosphi = cos(0.33), sinphi = sin(0.33);
00061 real  local_offset[3];
00062 float delay_time;
00063 float r_factor = 1.0;
00064 
00065 char  graph3d = 1, track = 0, cenergy = 0; // Convenient to distinguish these.
00066 
00067 char   temp_buffer[255];        // Convenient to share this temporary space.
00068 
00069 //-------------------------------------------------------------------------
00070 
00071 // base_point_size defines the size of the default "dot" size on the display.
00072 // It is used (with scaling) for stars and (without scaling) for nodes.
00073 
00074 local void set_base_point_size(float rel_point_size)
00075 {
00076     base_point_size = 2 * SMALL_DOT_SIZE;       // Default point size
00077     if (rel_point_size > 0.0)  base_point_size *= rel_point_size;
00078 }
00079 
00080 local float get_point_size(sdyn3* bi)
00081 {
00082     // Scaling by mass is too extreme.  Use sqrt(mass).
00083     // Also, impose a lower limit on the point size.
00084 
00085     if (point_scale_mode == 1)
00086         return max(SMALL_DOT_SIZE, base_point_size * sqrt(bi->get_mass()));
00087     else if (point_scale_mode == 2)
00088         return max(SMALL_DOT_SIZE, base_point_size * bi->get_radius());
00089     else
00090         return base_point_size;
00091 }
00092 
00093 local void draw_star_point(unsigned long win, float r, float s,
00094                            float actual_point_size, bool f_flag)
00095 
00096 // Plot a point of size actual_point_size at (r,s), filled or unfilled
00097 // according to the value of f_flag.
00098 
00099 {
00100     if (f_flag)
00101         lux_fill_arcf(win, r - actual_point_size/2, s - actual_point_size/2,
00102                       actual_point_size, actual_point_size, 0.0, 360.0);
00103     else
00104         lux_draw_arcf(win, r - actual_point_size/2, s - actual_point_size/2, 
00105                       actual_point_size, actual_point_size, 0.0, 360.0);
00106 }
00107 
00108 local void draw_links_2d(sdyn3* b, float r, float s)
00109 {
00110     for_all_daughters(sdyn3, b, bb) {
00111 
00112         float rr = bb->get_pos()[kx]-local_offset[kx];
00113         float ss = bb->get_pos()[ky]-local_offset[ky];
00114         float stemp;
00115 
00116         // Deal with the 9 possible locations of (rr,ss) with respect
00117         // to the box defined by xmin, xmax, ymin, and ymax.
00118 
00119         if (rr < xmin) {
00120 
00121             if (ss < ymin) {
00122                 float stemp = interp_to_x(r, s, rr, ss, xmin);
00123                 if (stemp >= ymin)
00124                     rr = xmin, ss = stemp;
00125                 else
00126                     rr = interp_to_y(r, s, rr, ss, ymin), ss = ymin;
00127             } else if (ss > ymax) {
00128                 if ((stemp = interp_to_x(r, s, rr, ss, xmin)) <= ymax)
00129                     rr = xmin, ss = stemp;
00130                 else
00131                     rr = interp_to_y(r, s, rr, ss, ymax), ss = ymax;
00132             } else
00133                 ss = interp_to_x(r, s, rr, ss, xmin), rr = xmin;
00134 
00135         } else if (rr > xmax) {
00136 
00137             if (ss < ymin) {
00138                 if ((stemp = interp_to_x(r, s, rr, ss, xmax)) >= ymin)
00139                     rr = xmax, ss = stemp;
00140                 else
00141                     rr = interp_to_y(r, s, rr, ss, ymin), ss = ymin;
00142             } else if (ss > ymax) {
00143                 if ((stemp = interp_to_x(r, s, rr, ss, xmax)) <= ymax)
00144                     rr = xmax, ss = stemp;
00145                 else
00146                     rr = interp_to_y(r, s, rr, ss, ymax), ss = ymax;
00147             } else
00148                 ss = interp_to_x(r, s, rr, ss, xmax), rr = xmax;
00149 
00150         } else {
00151 
00152             if (ss < ymin)
00153                 rr = interp_to_y(r, s, rr, ss, ymin), ss = ymin;
00154             else if (ss > ymax)
00155                 rr = interp_to_y(r, s, rr, ss, ymax), ss = ymax;
00156 
00157         }
00158         lux_draw_linef(win, r, s, rr, ss);
00159     }
00160 }                                          
00161 
00162 local void draw_links_3d(sdyn3* b, float r, float s)
00163 {
00164     for_all_daughters(sdyn3, b, bb) {
00165 
00166         float X = (float)bb->get_pos()[0] - local_offset[0] - origin[0];
00167         float Y = (float)bb->get_pos()[1] - local_offset[1] - origin[1];
00168         float Z = (float)bb->get_pos()[2] - local_offset[2] - origin[2];
00169 
00170         // Don't attempt to deal with the 27 possible locations
00171         // of (X, Y, Z) with respect to the box!
00172 
00173         float rr, ss;
00174         project3d(X, Y, Z, rr, ss, costheta, sintheta, cosphi, sinphi);
00175         lux_draw_linef(win, r, s, rr, ss);
00176     }
00177 }                                          
00178 
00179 local int plot_stars(sdyn3* b, int f_flag)
00180 {
00181     int  n_stars = b->n_leaves();
00182     if (b->get_oldest_daughter() == NULL) n_stars = 0;
00183 
00184     int  n_nodes = count_nodes(b);
00185     float r, s;
00186 
00187     if (graph3d) {
00188 
00189         for_all_nodes(sdyn3, b, bi) if (root || bi->get_parent()) {
00190             if (nodes || (bi->get_oldest_daughter() == NULL)) {
00191                 
00192                 float X = (float)bi->get_pos()[0] - local_offset[0] - origin[0];
00193                 float Y = (float)bi->get_pos()[1] - local_offset[1] - origin[1];
00194                 float Z = (float)bi->get_pos()[2] - local_offset[2] - origin[2];
00195                 float actual_point_size = get_point_size(bi);
00196 
00197 
00198                 if (   (X > (-lmax3d + actual_point_size))
00199                     && (X < ( lmax3d - actual_point_size)) 
00200                     && (Y > (-lmax3d + actual_point_size))
00201                     && (Y < ( lmax3d - actual_point_size)) 
00202                     && (Z > (-lmax3d + actual_point_size))
00203                     && (Z < ( lmax3d - actual_point_size))) {
00204 
00205                     // Should really sort by depth here...
00206 
00207                     project3d(X, Y, Z, r, s, costheta, sintheta,
00208                               cosphi, sinphi);
00209 
00210                     // Determine the color to use.
00211 
00212                     bool temp_flag = f_flag;
00213                     if (bi->get_oldest_daughter()) {
00214                         lux_set_color(win, lux_lookup_color(win, "grey"));
00215                         temp_flag =  1 - f_flag;
00216                     } else {
00217                         if (cenergy) {
00218 
00219                             char c;
00220 
00221                             compute_energies(b, bi, c);
00222                             lux_set_color(win, c_energy[c]);
00223 
00224                         } else if (bi->get_index() > 0) {
00225 
00226                             // Wrap the color map.
00227 
00228                             int ii = bi->get_index();
00229                             while (ii > N_COLORS) ii -= N_COLORS;
00230 
00231                             lux_set_color(win,c_index[ii]);
00232                         }
00233                     }
00234                     
00235                     if (track) 
00236                         lux_draw_pointf(win, r, s);
00237                     else
00238                         draw_star_point(win, r, s, actual_point_size, f_flag);
00239 
00240                     if (links && bi->get_oldest_daughter())
00241                         draw_links_3d(bi, r, s);
00242                 }
00243             }
00244         }
00245 
00246     } else {
00247 
00248         // Make a list of nodes, *including* the root node if root = 1.
00249         // Root will be at the start of the list if it is being displayed.
00250 
00251         sdyn3ptr* p = new sdyn3ptr[n_nodes+root];
00252 
00253         int ip = 0;
00254         for_all_nodes(sdyn3, b, bi) if (root || bi != b) p[ip++] = bi;
00255         if (ip != n_nodes+root) {
00256             cerr << "plot_stars: n_nodes = " << n_nodes+root
00257                 << " counted " << ip << endl;
00258             exit(1);
00259         }
00260 
00261         // Sort by kproj (note that kproj = 1, 2, or 3 for x, y, or z):
00262 
00263         for (ip = 0; ip < n_nodes+root; ip++)
00264             for (int jp = ip+1; jp < n_nodes+root; jp++)
00265                 if (p[jp]->get_pos()[kproj-1] < p[ip]->get_pos()[kproj-1]) {
00266                     sdyn3ptr bb = p[jp];
00267                     p[jp] = p[ip];
00268                     p[ip] = bb;
00269                 }
00270 
00271         // Plot ordered by depth.
00272 
00273         for (ip = 0; ip < n_nodes+root; ip++) {
00274             sdyn3 * bi = p[ip];
00275             if ( (root || bi != b)
00276                 && (nodes || bi->get_oldest_daughter() == NULL)) {
00277 
00278                 r = (float)bi->get_pos()[kx] - local_offset[kx];
00279                 s = (float)bi->get_pos()[ky] - local_offset[ky];
00280                 float actual_point_size = get_point_size(bi);
00281 
00282                 if (   (r > (xmin + actual_point_size))
00283                     && (r < (xmax - actual_point_size)) 
00284                     && (s > (ymin + actual_point_size))
00285                     && (s < (ymax - actual_point_size)) ) {
00286 
00287                     // Determine the color to use.
00288 
00289                     bool temp_flag = f_flag;
00290                     if (bi->get_oldest_daughter()) {
00291                         lux_set_color(win, lux_lookup_color(win, "grey"));
00292                         temp_flag =  1 - f_flag;
00293                     } else {
00294                         if (cenergy) {
00295 
00296                             char c;
00297 
00298                             compute_energies(b, bi, c);
00299                             lux_set_color(win, c_energy[c]);
00300 
00301                         } else if (bi->get_index() > 0) {
00302 
00303                             // Wrap the color map.
00304 
00305                             int ii = bi->get_index();
00306                             while (ii > N_COLORS) ii -= N_COLORS;
00307 
00308                             lux_set_color(win,c_index[ii]);
00309                         }
00310                     }
00311 
00312                     if (track) 
00313                         lux_draw_pointf(win, r, s);
00314                     else
00315                         draw_star_point(win, r, s,
00316                                         actual_point_size, temp_flag);
00317 
00318                     if (links && bi->get_oldest_daughter())
00319                         draw_links_2d(bi, r, s);
00320                 }
00321             }
00322         }
00323         delete p;
00324     }
00325     return n_stars;
00326 }
00327 
00328 local int type(int which)
00329 {
00330     if (which == tracking || which == graph3dim || which == colorenergy)
00331         return BUTTON_WINDOW;
00332     else
00333         return INPUT_WINDOW;
00334 }
00335 
00336 local int subtype(int which)
00337 {
00338     if (which == tracking || which == graph3dim || which == colorenergy)
00339         return CHECK_BUTTON;
00340     else
00341         return NO_TYPE;
00342 }
00343 
00344 local void set_temp_buffer(int which)
00345 {
00346     // Translate ID into a string representing the value.
00347 
00348     char c = -1;
00349     int i = -1;
00350     float f = VERY_LARGE_NUMBER;
00351 
00352     switch (which) {
00353         case colorenergy:       c = cenergy; break;
00354         case tracking:          c = track;   break;
00355         case graph3dim:         c = graph3d; break;
00356         case xminimum:          f = xmin; break;
00357         case xmaximum:          f = xmax; break;
00358         case yminimum:          f = ymin; break;
00359         case ymaximum:          f = ymax; break;
00360         case basepointsize:     f = base_point_size; break;
00361         case pointscalemode:    i = point_scale_mode; break;
00362         case lmax3D:            f = 2*lmax3d; break;
00363         case theta3D:           f = theta; break;
00364         case phi3D:             f = phi; break;
00365         case DelayTime:         f = delay_time; break;
00366         case dtheta3D:          f = dtheta; break;
00367         case Xorigin:           f = origin[0]; break;
00368         case Yorigin:           f = origin[0]; break;
00369         case Zorigin:           f = origin[0]; break;
00370         case view2D:            i = kproj; break;
00371         case originstar:        i = origin_star; break;
00372         default:                cerr << "xstarplot: diag error...\n";
00373     }
00374 
00375     if (f < VERY_LARGE_NUMBER)
00376         sprintf(temp_buffer, "%f", f);
00377     else if (c == -1)
00378         sprintf(temp_buffer, " %d ", i);
00379     else
00380         temp_buffer[0] = c;
00381 }
00382 
00383 local void set_diag_item(int which, char* id, int label_on_left,
00384                          int x1, int x2, int y)
00385 {
00386     set_temp_buffer(which);
00387 
00388     if (label_on_left) {        // label to left of (real) input box
00389 
00390         lux_set_item(dia, which, TEXT_WINDOW, NO_TYPE,
00391                      _R_(x1), _R_(y), strlen(id), id);
00392         lux_set_item(dia, which, INPUT_WINDOW, NO_TYPE,
00393                      _R_(x2), _R_(y), 10, temp_buffer);
00394 
00395     } else {                    // label to right of (integer) input box/button
00396 
00397         if (type(which) == INPUT_WINDOW) {
00398 
00399             lux_set_item(dia, which, INPUT_WINDOW, NO_TYPE,
00400                          _R_(x1), _R_(y), 3, temp_buffer);
00401             lux_set_item(dia, which, TEXT_WINDOW, NO_TYPE,
00402                          _R_(x2), _R_(y), strlen(id), id);
00403 
00404         } else {
00405 
00406             lux_set_item(dia, which, BUTTON_WINDOW, CHECK_BUTTON,
00407                          _R_(x1), _R_(y), 1, temp_buffer);
00408             lux_set_item(dia, which, TEXT_WINDOW, CHECK_BUTTON,
00409                          _R_(x2), _R_(y), strlen(id), id);
00410         }
00411     }
00412 }
00413 
00414 local void initialize_dialog(int xorigin, int yorigin)
00415 {
00416     // Initialize dialog box material (note: y is measured up from bottom!)
00417 
00418     float save = r_factor;
00419     if (r_factor > 1) r_factor = 1;
00420 
00421     int xsize = 550;
00422     int ysize = 550;
00423     if (r_factor <= 0.6) xsize = (int) (xsize / (r_factor/0.6));
00424     dia = lux_open_dialog(xorigin, yorigin, _R_(xsize), _R_(ysize));
00425 
00426     // ---------- 3-D origin info across top of box (special case): ----------
00427 
00428      lux_set_item(dia, Origin, TEXT_WINDOW, NO_TYPE,
00429                  _R_(70), _R_(500), 6, "origin");
00430 
00431     sprintf(temp_buffer, "%f", origin[0]);
00432     lux_set_item(dia, Xorigin, INPUT_WINDOW, NO_TYPE,
00433                  _R_(160), _R_(500), 10, temp_buffer);
00434 
00435     sprintf(temp_buffer, "%f", origin[1]);
00436     lux_set_item(dia, Yorigin, INPUT_WINDOW, NO_TYPE,
00437                  _R_(280), _R_(500), 10, temp_buffer);
00438 
00439     sprintf(temp_buffer, "%f", origin[2]);
00440     lux_set_item(dia, Zorigin, INPUT_WINDOW, NO_TYPE,
00441                  _R_(400), _R_(500), 10, temp_buffer);
00442 
00443 
00444     // ---------- Left-hand column of dialog box: ----------
00445 
00446     //          xmin
00447     //          xmax
00448     //          ymin
00449     //          ymax
00450     //          3-D cube size
00451     //          projection theta
00452     //          projection phi
00453     //          projection dtheta
00454     //          point scale
00455     //          delay time
00456 
00457     // Minima and maxima of 2-D display:
00458 
00459     set_diag_item(xminimum, "xmin", 1, 70, 160, 460);
00460     set_diag_item(xmaximum, "xmax", 1, 70, 160, 420);
00461     set_diag_item(yminimum, "ymin", 1, 70, 160, 380);
00462     set_diag_item(ymaximum, "ymax", 1, 70, 160, 340);
00463 
00464     // 3-D cube size:
00465 
00466     set_diag_item(lmax3D, "cube size", 1, 70, 180, 300);
00467 
00468     // 3-D projection direction:
00469 
00470     set_diag_item(theta3D,  "theta",  1, 70, 180, 260);
00471     set_diag_item(phi3D,    "phi",    1, 70, 180, 220);
00472     set_diag_item(dtheta3D, "dtheta", 1, 70, 180, 180);
00473 
00474     // Point size:
00475 
00476     set_diag_item(basepointsize, "point scale", 1, 70, 180, 140);
00477 
00478     // Delay time:
00479 
00480     set_diag_item(DelayTime, "delay time", 1, 70, 180, 100);
00481 
00482 
00483     // ---------- Right-hand column of dialog box: ----------
00484 
00485     //          origin star
00486     //          2-D view axis
00487     //          3-D plot?
00488 
00489     //          point display mode
00490     //          color by energy?
00491     //          track particles?
00492 
00493     // Origin star:
00494 
00495     set_diag_item(originstar, "Origin Star", 0, 320, 360, 460);
00496 
00497     // 2-D view axis:
00498 
00499     set_diag_item(view2D, "2-D view axis", 0, 320, 360, 420);
00500 
00501     // 3D graphics:
00502 
00503     set_diag_item(graph3dim, "3-D graph", 0, 320, 350, 380);
00504 
00505     // Point display mode
00506 
00507     set_diag_item(pointscalemode, "point display mode", 0, 320, 360, 300);
00508 
00509     // Color by energy:
00510 
00511     set_diag_item(colorenergy, "color by energy", 0, 320, 350, 260);
00512 
00513     // Track particles:
00514 
00515     set_diag_item(tracking, "track particles", 0, 320, 350, 220);
00516 
00517     // lux_draw_palette(dia);
00518 
00519     // OK/cancel:
00520 
00521     lux_set_item(dia, ok, BUTTON_WINDOW, OK_BUTTON,
00522                  _R_(100), _R_(25), 9, "OK, CLEAR");
00523     lux_set_item(dia, ok_keep, BUTTON_WINDOW, OK_KEEP_BUTTON,
00524                  _R_(250), _R_(25), 8, "OK, KEEP");
00525     lux_set_item(dia, cancel, BUTTON_WINDOW, CANCEL_BUTTON,
00526                  _R_(400), _R_(25), 6, "CANCEL");
00527 
00528     r_factor = save;
00529 }
00530 
00531 local void make_relative_to_root(sdyn3* b)
00532 {
00533     vector root_pos = b->get_pos();
00534     vector root_vel = b->get_vel();
00535     for_all_nodes(sdyn3, b, bi) {
00536         bi->inc_pos(-root_pos);
00537         bi->inc_vel(-root_vel);
00538     }
00539 }
00540 
00541 local void update_diag_item(int which)
00542 {
00543     set_temp_buffer(which);
00544     lux_update_itemvalue(dia, which, type(which), subtype(which),
00545                          temp_buffer);
00546 }
00547 
00548 local void get_diag_string(int which) {
00549 }
00550 
00551 // Overloaded function:
00552 
00553 local void read_diag_item(int which, float& f)
00554 {
00555     lux_get_itemvalue(dia, which, type(which), subtype(which),
00556                       temp_buffer);
00557     sscanf(temp_buffer, "%f", &f);
00558 }
00559 
00560 local void read_diag_item(int which, int& i)
00561 {
00562     lux_get_itemvalue(dia, which, type(which), subtype(which),
00563                       temp_buffer);
00564     sscanf(temp_buffer, "%d", &i);
00565 }
00566 
00567 local void read_diag_item(int which, char& c)
00568 {
00569     lux_get_itemvalue(dia, which, type(which), subtype(which),
00570                       temp_buffer);
00571     c = temp_buffer[0];
00572 }
00573 
00574 local void update_from_dialog(bool b_flag)
00575 {
00576     // Read all data from the dialog box.  Note that it is
00577     // important that the box be maintained correctly.
00578     // Any errors in the data displayed in the box will be
00579     // propogated to the "real" data by this procedure.
00580 
00581     read_diag_item(view2D, kproj);
00582     if (kproj == 1)
00583         kx = 1, ky = 2;
00584     else if (kproj == 2)
00585         kx = 2, ky = 0;
00586     else
00587         kx = 0, ky = 1;
00588 
00589     read_diag_item(originstar, origin_star);
00590 
00591     read_diag_item(Xorigin, origin[0]);
00592     read_diag_item(Yorigin, origin[1]);
00593     read_diag_item(Zorigin, origin[2]);
00594 
00595     read_diag_item(lmax3D, lmax3d);
00596     lmax3d /= 2;                        // Display twice lmax3d
00597 
00598     // ------------------------------------------------------
00599 
00600     // NOTE: xmin, xmax, ymin, ymax are IRRELEVANT, since they
00601     // will be forced to be consistent with origin and lmax3d.
00602 
00603     read_diag_item(xminimum, xmin);
00604     read_diag_item(xmaximum, xmax);
00605     read_diag_item(yminimum, ymin);
00606     read_diag_item(ymaximum, ymax);
00607 
00608     set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
00609 
00610     // Make sure entries in the dialog box are consistent:
00611 
00612     update_diag_item(xminimum);
00613     update_diag_item(xmaximum);
00614     update_diag_item(yminimum);
00615     update_diag_item(ymaximum);
00616 
00617     // ------------------------------------------------------
00618 
00619     read_diag_item(basepointsize,base_point_size);
00620 
00621     read_diag_item(theta3D, theta);
00622     costheta = cos(theta);
00623     sintheta = sin(theta);
00624 
00625     read_diag_item(phi3D, phi);
00626     cosphi = cos(phi);
00627     sinphi = sin(phi);
00628 
00629     read_diag_item(dtheta3D, dtheta);
00630 
00631     read_diag_item(DelayTime, delay_time);
00632 
00633     read_diag_item(colorenergy, cenergy);
00634     show_color_scheme(colwin, c_energy, c_index,
00635                       r_factor, cenergy, b_flag, 1);
00636 
00637     read_diag_item(tracking, track);
00638     read_diag_item(graph3dim, graph3d);
00639 
00640     lux_clear_window(win);
00641 
00642     if (graph3d) {
00643         lux_setup_axis(win, -FAC3D*lmax3d, FAC3D*lmax3d, 
00644                             -FAC3D*lmax3d, FAC3D*lmax3d);
00645         draw3d_axis(win, lmax3d, costheta, sintheta,
00646                     cosphi, sinphi);
00647     } else {
00648         lux_setup_axis(win, xmin, xmax, ymin, ymax);
00649         draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
00650     }
00651 
00652     read_diag_item(pointscalemode, point_scale_mode);
00653 }
00654 
00655 local void show_static_rotation(sdyn3* b, bool f_flag)
00656 {
00657     // Handle static rotation using lux_check_keypress:
00658 
00659     show_instructions(instr, r_factor,
00660                       "Pause...\n  r: rotate, (p,c): continue                ",
00661                       1);
00662 
00663     char rotate = 0;
00664     do {
00665         if (lux_check_keypress(win,'r')) {
00666             show_instructions(instr, r_factor,
00667                       "Rotate... p: pause, c: continue, c: exit rotation     ",
00668                               0);
00669             rotate = 1;
00670 
00671             do {
00672                 theta += dtheta;
00673                 if (theta < 0.0)
00674                     theta += 2*PI;
00675                 else if (theta > PI*2)
00676                     theta -= 2*PI;
00677                 costheta = cos(theta);
00678                 sintheta = sin(theta);
00679                 update_diag_item(theta3D);
00680 
00681                 lux_clear_current_region(win);
00682                 draw3d_axis(win, lmax3d, costheta, sintheta,
00683                             cosphi, sinphi);
00684 
00685                 for_all_leaves(sdyn3, b, bi) {
00686                     float X, Y, Z;
00687                     X = (float)bi->get_pos()[0] - origin[0];
00688                     Y = (float)bi->get_pos()[1] - origin[1];
00689                     Z = (float)bi->get_pos()[2] - origin[2];
00690                     float actual_point_size
00691                         = get_point_size(bi);
00692 
00693                     if (   (X > (-lmax3d + actual_point_size))
00694                         && (X < ( lmax3d - actual_point_size)) 
00695                         && (Y > (-lmax3d + actual_point_size))
00696                         && (Y < ( lmax3d - actual_point_size)) 
00697                         && (Z > (-lmax3d + actual_point_size))
00698                         && (Z < ( lmax3d - actual_point_size))){
00699 
00700                         float r, s;
00701                         project3d(X, Y, Z, r, s,
00702                                   costheta, sintheta,
00703                                   cosphi, sinphi);
00704 
00705                         if (cenergy) {
00706                             char c;
00707                             compute_energies(b, bi, c);
00708                             lux_set_color(win,c_energy[c]);
00709                         } else if (bi->get_index()>0) {
00710                             if (bi->get_index() <= N_COLORS)
00711                                 lux_set_color(win,
00712                                               c_index[bi->get_index()]);
00713                             else
00714                                 lux_set_color(win,c_index[1]);
00715                         }
00716 
00717                         if (track) 
00718                             lux_draw_pointf(win,r,s);
00719                         else
00720                             if (f_flag)
00721                                 lux_fill_arcf(win,
00722                                               r - actual_point_size/2,
00723                                               s - actual_point_size/2, 
00724                                               actual_point_size,
00725                                               actual_point_size,
00726                                               0.0, 360.0);
00727                             else
00728                                 lux_draw_arcf(win,
00729                                               r - actual_point_size/2,
00730                                               s - actual_point_size/2, 
00731                                               actual_point_size,
00732                                               actual_point_size,
00733                                               0.0, 360.0);
00734                     }
00735                 }
00736                 update_with_delay(win, delay_time);
00737 
00738                 if (lux_check_keypress(win,'p')) 
00739                     while(!lux_check_keypress(win,'c'));
00740 
00741             } while(!lux_check_keypress(win,'c'));
00742 
00743             // Flush "c"s from input stream:
00744 
00745             while(lux_check_keypress(win,'c'));
00746         }
00747 
00748     } while(!lux_check_keypress(win,'p')
00749             && !lux_check_keypress(win,'c') && !rotate);
00750 
00751     rotate = 0;
00752     lux_set_color(win,c_energy[default_color]);
00753 }
00754 
00755 local void check_for_input(unsigned long win, sdyn3* b,
00756                            bool b_flag, bool f_flag)
00757 {
00758     // Check for interactive input.
00759     // Loop through the input buffer until no more events remain.
00760 
00761     // lux_next_keypress gets and discards the next key in the input stream.
00762 
00763     char key, string[20], shift, control;
00764 
00765     while (lux_next_keypress(win, &key, string, &shift, &control)) {
00766 
00767         // A "defined" key has been pressed.  See if it is one we want.
00768 
00769         if (key == 0                    // key = 0 for non-ASCII character
00770                                         // e.g. Up, Down, Home, etc.--see win.c
00771             || key == 'h'
00772             || key == 'a') {
00773 
00774             // "Shift" functions:
00775 
00776             int change = 0;
00777 
00778             if (key == 'h')             // Keyboard lookalikes for
00779                 change = 4;             // function keys.
00780             else if (key == 'a')
00781                 change = 5;
00782 
00783             else if (strcmp(string, "Up") == 0) {
00784                 change = ky + 1;
00785                 origin[ky] += lmax3d/2;
00786             } else if (strcmp(string, "Down") == 0) {
00787                 change = ky + 1;
00788                 origin[ky] -= lmax3d/2;
00789             } else if (strcmp(string, "Right") == 0) {
00790                 change = kx + 1;
00791                 origin[kx] += lmax3d/2;
00792             } else if (strcmp(string, "Left") == 0) {
00793                 change = kx + 1;
00794                 origin[kx] -= lmax3d/2;
00795             } else if (strcmp(string, "PgUp") == 0) {
00796                 change = kproj;
00797                 origin[kproj-1] += lmax3d/2;
00798             } else if (strcmp(string, "PgDn") == 0) {
00799                 change = kproj;
00800                 origin[kproj-1] -= lmax3d/2;
00801             } else if (strcmp(string, "Home") == 0)
00802                 change = 4;
00803             else if (strcmp(string, "R11") == 0)
00804                 change = 5;
00805 
00806             if (change) {
00807 
00808                 if (change == 4) {
00809                     origin[kx] = origin[ky] = origin[kproj-1] = 0;
00810                 } else if (change == 5) {
00811 
00812                     real save = lmax3d;
00813                     lmax3d = 0;
00814 
00815                     // Use all dimensions in redetermining the scale!
00816 
00817                     for_all_leaves(sdyn3, b, bi)
00818                         for (int kk = 0; kk < 3; kk++)
00819                             lmax3d = max(lmax3d, abs(bi->get_pos()[kk]
00820                                                      - local_offset[kk]));
00821                 
00822                     // Round lmax3d up to something reasonable:
00823                 
00824                     lmax3d *= 1.2;
00825                     if (lmax3d <= 0) lmax3d = 1;
00826                 
00827                     real m_log = log10(lmax3d);
00828                     int i_log = (int) m_log;
00829                     if (m_log - i_log > 0.699) i_log++;
00830                     real scale = pow(10.0, i_log);
00831                     lmax3d = ((int) (lmax3d/scale) + 1) * scale;
00832 
00833                     origin[kx] = origin[ky] = origin[kproj-1] = 0;              
00834                     base_point_size *= lmax3d/save;
00835                 }
00836 
00837                 set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
00838 
00839                 // Update all entries in the dialog box, as needed:
00840 
00841                 update_diag_item(xminimum);
00842                 update_diag_item(xmaximum);
00843                 update_diag_item(yminimum);
00844                 update_diag_item(ymaximum);
00845 
00846                 if (change == 1 || change >= 4) update_diag_item(Xorigin);
00847                 if (change == 2 || change >= 4) update_diag_item(Yorigin);
00848                 if (change == 3 || change >= 4) update_diag_item(Zorigin);
00849 
00850                 if (change == 5) {
00851                     update_diag_item(basepointsize);
00852                     update_diag_item(lmax3D);
00853                 }
00854 
00855                 // No redrawing of axes is needed for 3D graphs for change != 5.
00856                 
00857                 if (!graph3d) {
00858                     lux_clear_window(win);
00859                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
00860                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
00861                 } else if (change == 5) {
00862                     lux_clear_window(win);
00863                     lux_setup_axis(win, -FAC3D*lmax3d, FAC3D*lmax3d, 
00864                                    -FAC3D*lmax3d, FAC3D*lmax3d);
00865                     draw3d_axis(win, lmax3d, costheta, sintheta,
00866                                 cosphi, sinphi);
00867                 }
00868             }
00869 
00870         } else {                        // Key is the ASCII code.
00871 
00872             // Check for "quit" first.
00873 
00874             if (key == 'q') {
00875                 show_instructions(instr, r_factor,
00876                     "\n   Quitting...                                        ",
00877                                   1);
00878                 lux_exit();     // Note that this is now quick_exit
00879             }
00880 
00881             //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
00882 
00883             if (graph3d && !track) {
00884                 int button;
00885 
00886                 // 3-D operations:
00887 
00888                 if (key== 'R') {
00889 
00890                     show_static_rotation(b, f_flag);
00891 
00892                 } else {
00893 
00894                     button = lux_check_buttonpress(win);
00895 
00896                     if (key == '^') {
00897                         phi = phi - dtheta;
00898                         if (phi < -PI) phi = phi + 2.0*PI;
00899                         cosphi = cos(phi);
00900                         sinphi = sin(phi);
00901                         update_diag_item(phi3D);
00902                     } else if (key == 'V') {
00903                         phi = phi + dtheta;
00904                         if (phi > PI) phi = phi - 2.0*PI;
00905                         cosphi = cos(phi);
00906                         sinphi = sin(phi);
00907                         update_diag_item(phi3D);
00908                     } else if (key == '<') {
00909                         theta = theta + dtheta;
00910                         if (theta > 2*PI) theta -= 2*PI;
00911                         costheta = cos(theta);
00912                         sintheta = sin(theta);
00913                         update_diag_item(theta3D);
00914                     } else if (key == '>') {
00915                         theta = theta - dtheta;
00916                         if (theta < 0.0) theta += 2*PI;
00917                         costheta = cos(theta);
00918                         sintheta = sin(theta);
00919                         update_diag_item(theta3D);
00920                     }
00921                 }
00922             }
00923 
00924             //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
00925 
00926             if (key == 'p' || key == 'P') {     // Modify base point size
00927 
00928                 if (key == 'p') {
00929                     if (base_point_size > 4*lmax3d/win_size)
00930                         base_point_size /= PFAC;
00931                 } else
00932                     base_point_size *= PFAC;
00933 
00934                 update_diag_item(basepointsize);
00935 
00936             } else if (key == 'n') {            // Toggle tree display
00937 
00938                 nodes = 1 - nodes;
00939                 if (nodes == 0) links = 0;
00940 
00941             } else if (key == 'l') {            // Toggle links/nodes
00942 
00943                 links = 1 - links;
00944                 if (links == 1) nodes = 1;
00945 
00946             } else if (key == 'r') {            // Show root node
00947 
00948                 root = 1 - root;
00949 
00950             } else if (key == 't') {            // Toggle tracking:
00951 
00952                 track = 1 - track;
00953                 update_diag_item(tracking);
00954 
00955             } else if (key == '3' && !graph3d) {            // Enter 3D mode
00956 
00957                 graph3d = 1;
00958                 update_diag_item(graph3dim);
00959                 lux_clear_window(win);
00960                 lux_setup_axis(win, -FAC3D*lmax3d, FAC3D*lmax3d, 
00961                                -FAC3D*lmax3d, FAC3D*lmax3d);
00962                 draw3d_axis(win, lmax3d, costheta, sintheta, cosphi, sinphi);
00963 
00964             } else if (graph3d &&
00965                        (key == '2' || key == 'x'
00966                           || key == 'y' || key == 'k')) {   // Enter 2D mode
00967 
00968                 graph3d = 0;
00969                 update_diag_item(graph3dim);
00970                 lux_clear_window(win);
00971                 lux_setup_axis(win, xmin, xmax, ymin, ymax);
00972                 draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
00973 
00974             } else if (key == 'e') {            // Color by energy
00975 
00976                 cenergy = !cenergy;
00977                 update_diag_item(colorenergy);
00978                 show_color_scheme(colwin, c_energy, c_index,
00979                                   r_factor, cenergy, b_flag, 1);
00980 
00981             } else if (key == '+' || key == '=') {      // Increase delay time
00982 
00983                 delay_time += 0.05;
00984                 update_diag_item(DelayTime);
00985 
00986             } else if (key == '-') {            // Decrease delay time
00987 
00988                 delay_time -= 0.05;
00989                 if (delay_time < 0.0) delay_time = 0.0;
00990                 update_diag_item(DelayTime);
00991 
00992             } if (key == '0') {                 // Set delay time = 0
00993 
00994                 delay_time = 0.0;
00995                 update_diag_item(DelayTime);
00996 
00997                 // Flush the queue:
00998 
00999                 while(lux_check_keypress(win,'+')
01000                       || lux_check_keypress(win,'-'));
01001 
01002             } else if (key == 'z') {            // Zoom in
01003                                                 // (fixed axes and origin)
01004 
01005                 lmax3d /= ZOOM;
01006                 set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01007                 base_point_size /= ZOOM;
01008 
01009                 if (graph3d) {
01010                     lux_clear_window(win);
01011                     lux_setup_axis(win,-FAC3D*lmax3d,FAC3D*lmax3d, 
01012                                    -FAC3D*lmax3d, FAC3D*lmax3d);
01013                     draw3d_axis(win, lmax3d, costheta, sintheta,
01014                                 cosphi, sinphi);
01015                 } else {
01016                     lux_clear_window(win);
01017                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
01018                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01019                 }
01020 
01021                 update_diag_item(basepointsize);
01022                 update_diag_item(lmax3D);
01023                 update_diag_item(xminimum);
01024                 update_diag_item(xmaximum);
01025                 update_diag_item(yminimum);
01026                 update_diag_item(ymaximum);
01027 
01028             } else if (key == 'Z') {            // Zoom out
01029 
01030                 lmax3d *= ZOOM;
01031                 set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01032                 base_point_size *= ZOOM;
01033 
01034                 if (graph3d) {
01035                     lux_clear_window(win);
01036                     lux_setup_axis(win, -FAC3D*lmax3d, FAC3D*lmax3d, 
01037                                    -FAC3D*lmax3d, FAC3D*lmax3d);
01038                     draw3d_axis(win, lmax3d, costheta, sintheta,
01039                                 cosphi, sinphi);
01040                 } else {
01041                     lux_clear_window(win);
01042                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
01043                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01044                 }
01045 
01046                 update_diag_item(basepointsize);
01047                 update_diag_item(lmax3D);
01048                 update_diag_item(xminimum);
01049                 update_diag_item(xmaximum);
01050                 update_diag_item(yminimum);
01051                 update_diag_item(ymaximum);
01052 
01053             } else if (key == 'k') {            // Change projection axis to z
01054                                                 // (fixed origin, scale)
01055                 if (graph3d) {
01056 
01057                     theta = -PI/2;
01058                     phi = PI/2.0;
01059                     costheta = cos(theta);
01060                     sintheta = sin(theta);
01061                     cosphi = cos(phi);
01062                     sinphi = sin(phi);
01063 
01064                     update_diag_item(theta3D);
01065                     update_diag_item(phi3D);
01066 
01067                     while (lux_check_keypress(win,'<')
01068                            || lux_check_keypress(win,'>')
01069                            || lux_check_keypress(win,'^')
01070                            || lux_check_keypress(win,'V') );
01071 
01072                 } else {
01073 
01074                     kproj = 3;
01075                     kx = 0;
01076                     ky = 1;
01077                     set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01078                     lux_clear_window(win);
01079                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
01080                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01081 
01082                     update_diag_item(view2D);
01083 
01084                 }
01085 
01086             } else if (key == 'y') {            // Change projection axis to y
01087 
01088                 if (graph3d) {
01089 
01090                     theta = -PI/2.0;
01091                     phi = 0.0;
01092                     costheta = cos(theta);
01093                     sintheta = sin(theta);
01094                     cosphi = cos(phi);
01095                     sinphi = sin(phi);
01096 
01097                     while (lux_check_keypress(win,'<')
01098                            || lux_check_keypress(win,'>')
01099                            || lux_check_keypress(win,'^')
01100                            || lux_check_keypress(win,'V') );
01101 
01102                     update_diag_item(theta3D);
01103                     update_diag_item(phi3D);
01104 
01105                 } else {
01106 
01107                     kproj = 2;
01108                     kx = 2;
01109                     ky = 0;
01110                     set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01111                     lux_clear_window(win);
01112                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
01113                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01114 
01115                     update_diag_item(view2D);
01116 
01117                 }
01118 
01119             } else if (key == 'x') {            // Change projection axis to x
01120 
01121                 if (graph3d) {
01122                     theta = 0.0;
01123                     phi = 0.0;
01124                     costheta = cos(theta);
01125                     sintheta = sin(theta);
01126                     cosphi = cos(phi);
01127                     sinphi = sin(phi);
01128 
01129                     while (lux_check_keypress(win,'<')
01130                            || lux_check_keypress(win,'>')
01131                            || lux_check_keypress(win,'^')
01132                            || lux_check_keypress(win,'V') );
01133 
01134                     update_diag_item(theta3D);
01135                     update_diag_item(phi3D);
01136 
01137                 } else {
01138 
01139                     kproj = 1;
01140                     kx = 1;
01141                     ky = 2;
01142                     set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01143                     lux_clear_window(win);
01144                     lux_setup_axis(win, xmin, xmax, ymin, ymax);
01145                     draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01146 
01147                     update_diag_item(view2D);
01148 
01149                 }
01150             }
01151 
01152             // Select new origin: 'o' ==> point in space, 'O' ==> star
01153 
01154             if (key == 'o' && !graph3d) {
01155 
01156                 show_instructions(instr, r_factor,
01157                 " Use mouse-right to select new origin...\n",1);
01158                 float r, s;
01159                 get_mouse_position(win, &r, &s);
01160 
01161                 // Looks like r and s come back in the correct units!
01162 
01163                 origin[kx] = r;
01164                 origin[ky] = s;
01165                 set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01166 
01167                 lux_clear_window(win);
01168                 lux_setup_axis(win, xmin, xmax, ymin, ymax);
01169                 draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01170 
01171                 origin_star = -1;
01172 
01173                 update_diag_item(Xorigin);
01174                 update_diag_item(Yorigin);
01175                 update_diag_item(Zorigin);
01176                 update_diag_item(xminimum);
01177                 update_diag_item(xmaximum);
01178                 update_diag_item(yminimum);
01179                 update_diag_item(ymaximum);
01180                 update_diag_item(originstar);
01181 
01182             } else if (!graph3d && key == 'O') {
01183 
01184                 show_instructions(instr, r_factor,
01185                 " Use mouse-right to select new origin star...\n",1);
01186                 float r, s;
01187                 get_mouse_position(win, &r, &s);
01188 
01189                 // Looks like r and s come back in the correct units!
01190 
01191                 // Find the star closest to the (2-D) mouse position.
01192 
01193                 origin_star = nearest_index(b, r, s, kx, ky);
01194 
01195                 if (origin_star > 0) {
01196                     for (int kk = 0; kk < 3; kk++) origin[kk] = 0;
01197                     set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01198                 }
01199 
01200                 lux_clear_window(win);
01201                 lux_setup_axis(win, xmin, xmax, ymin, ymax);
01202                 draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01203 
01204                 update_diag_item(Xorigin);
01205                 update_diag_item(Yorigin);
01206                 update_diag_item(Zorigin);
01207                 update_diag_item(xminimum);
01208                 update_diag_item(xmaximum);
01209                 update_diag_item(yminimum);
01210                 update_diag_item(ymaximum);
01211                 update_diag_item(originstar);
01212             }
01213 
01214             // -------------------------------------------------------------
01215 
01216             while(lux_check_keypress(win,'r')); // Clean up any replay garbage
01217 
01218             // Idle mode and dialog box input:
01219 
01220             int return_value = 0;
01221 
01222             if (key == 'i') {
01223 
01224                 show_instructions(instr, r_factor,
01225                   " Idle...\n  d: dialog, r: replay, c: continue: q-quit",1);
01226                 return_value = lux_getevent();
01227 
01228             } else if (key == 'd') {
01229 
01230                 show_instructions(instr, r_factor,
01231 "Dialog mode keyboard options\n\n\
01232   Dialog window:\n\
01233     o   OK, keep dialog window\n\
01234     O   OK, close dialog window\n\
01235     c   CANCEL, keep dialog window\n\
01236     C   CANCEL, close dialog window\n\n\
01237   Other windows:\n\
01238     r   replay
01239     c   continue (keep dialog window,\n\
01240                   ignore dialog changes)\n\
01241     q   quit\n\
01242 ",1);
01243                 lux_show_dialog(dia);
01244                 return_value = lux_getevent();
01245             }
01246             
01247             // Update all values if OK button was clicked (clicking
01248             // OK forces a return value of 3).
01249             
01250             if (return_value == 3) update_from_dialog(b_flag);
01251         }
01252     }
01253 }
01254 
01255 /*-----------------------------------------------------------------------------
01256  *  xstarplot.C  --  project the positions of all particles onto the screen
01257  *                   input:   pn: a pointer to a nbody system,
01258  *                              k: the number of the coordinate axis along which
01259  *                                 the projection is directed, with {k,x,y}
01260  *                                 having a right-handed orientation,
01261  *                           xmax: displayed x-axis spans [-xmax, xmax]
01262  *                           ymax: displayed y-axis spans [-ymax, ymax]
01263  *-----------------------------------------------------------------------------
01264  */
01265 
01266 void xstarplot(sdyn3* b, float scale, int k, int d, float lmax,
01267                int point_mode, float rel_point_size,
01268                float D, int ce,
01269                bool b_flag, bool f_flag, bool t_flag,
01270                int init_flag)
01271 {
01272     r_factor = scale;
01273 
01274     if (init_flag == 0) {
01275 
01276         int xorigin, yorigin;
01277 
01278         if (t_flag) nodes = links = 1;
01279 
01280         // Open windows for plotting and instructions, set up color
01281         // schemes, etc.
01282 
01283         if (init_status == 0) initialize_graphics(r_factor, b_flag,
01284                                                   win, instr, colwin,
01285                                                   c_index, c_energy,
01286                                                   win_size, xorigin, yorigin);
01287         lux_clear_window(win);
01288         lux_update_fg(win);
01289         lux_clear_window(colwin);
01290         lux_update_fg(colwin);
01291         lux_clear_window(instr);
01292         lux_update_fg(instr);
01293 
01294         // Determine which axes to plot:
01295 
01296         kproj = k;
01297         switch (kproj) {
01298             case 1:     kx = 1; ky = 2; graph3d = 0; break;
01299             case 2:     kx = 2; ky = 0; graph3d = 0; break;
01300             case 3:     kx = 0; ky = 1; graph3d = 0; break;
01301             default:    cerr << "xstarplot: k = " << k
01302                              << ": illegal value; choose from {1, 2, 3}"
01303                              << endl;
01304                         exit(0);
01305         }
01306 
01307         switch(d) {
01308             case 2:     graph3d = 0;  break;
01309             case 3:     graph3d = 1;  break;
01310             default:    cerr << "xstarplot: d = " << d
01311                              << " illegal value; choose from {2, 3)"
01312                              << endl;
01313                         exit(0);
01314         }
01315 
01316         point_scale_mode = point_mode;
01317 
01318         cenergy = ce?1:0;
01319         show_color_scheme(colwin, c_energy, c_index,
01320                           r_factor, cenergy, b_flag,1);
01321 
01322         delay_time = D;
01323 
01324         if (lmax > 0.0)
01325 
01326             lmax3d = lmax;
01327 
01328         else {
01329 
01330             lmax3d = 0;
01331 
01332             // Use all dimensions in determining the initial scale!
01333 
01334             for_all_leaves(sdyn3, b, bi)
01335                 for (int kk = 0; kk < 3; kk++)
01336                     lmax3d = max(lmax3d, abs(bi->get_pos()[kk]));
01337 
01338             // Round lmax3d up to something reasonable:
01339 
01340             lmax3d *= 1.2;
01341             if (lmax3d <= 0) lmax3d = 1;
01342 
01343             real m_log = log10(lmax3d);
01344             int i_log = (int) m_log;
01345             if (m_log - i_log > 0.699) i_log++;
01346             real scale = pow(10.0, i_log);
01347             lmax3d = ((int) (lmax3d/scale) + 1) * scale;
01348         }
01349 
01350         for (int ko = 0; ko < 3; ko++) origin[ko] = 0;
01351 
01352         set_limits(origin, lmax3d, kx, xmin, xmax, ky, ymin, ymax);
01353         set_base_point_size(rel_point_size);
01354         
01355         if (graph3d) {
01356             lux_setup_axis(win, -FAC3D*lmax3d, FAC3D*lmax3d, 
01357                            -FAC3D*lmax3d, FAC3D*lmax3d);      
01358             draw3d_axis(win, lmax3d, costheta, sintheta, cosphi, sinphi);
01359         } else {
01360             lux_setup_axis(win, xmin, xmax, ymin, ymax);
01361             draw2d_axis(win, xmin, xmax, ymin, ymax, kproj);
01362         }
01363 
01364         // Complete the initialization with the dialog box.
01365 
01366         if (init_status == 0) initialize_dialog(xorigin, yorigin);
01367 
01368         init_status = 1;
01369             
01370     } else if (init_flag < 0) {
01371 
01372         show_instructions(instr, r_factor,
01373                   "End of Data.  Idle now.\n  r: replay, (c,q): quit", 1);
01374 
01375         lux_getevent();
01376         exit(0);
01377 
01378     }
01379 
01380     // lux_getevent();
01381     // exit(1);
01382 
01383     //------------------------------------------------------------------------
01384 
01385     // End of initialization.  Start by (re)drawing axes and instructions.
01386 
01387     if (track == 0) {
01388 
01389         lux_clear_current_region(win);
01390         if (graph3d) draw3d_axis(win, lmax3d, costheta, sintheta,
01391                                  cosphi, sinphi);
01392     }
01393 
01394     show_main_instructions(instr, r_factor, graph3d, 1);
01395 
01396     //------------------------------------------------------------------------
01397 
01398     // Always express all positions and velocities relative to the root node.
01399     // (There is no requirement that the root node be at rest at the origin...)
01400 
01401     make_relative_to_root(b);
01402 
01403     // Determine local offset, if any.
01404 
01405     if (origin_star > 0) {
01406         for_all_leaves(sdyn3, b, bi)
01407             if (bi->get_index() == origin_star) {
01408                 for (int kk = 0; kk < 3; kk++)
01409                     local_offset[kk] = bi->get_pos()[kk];
01410                 break;
01411             }
01412     } else
01413         for (int kk = 0; kk < 3; kk++) local_offset[kk] = 0;
01414 
01415     // Plot the data points (sorted by kproj-component in the 2D case):
01416 
01417     int n_stars = plot_stars(b, f_flag);
01418 
01419     // Update headers.
01420 
01421     lux_set_color(win, c_energy[default_color]);
01422 
01423     if (graph3d) {
01424 
01425         sprintf(temp_buffer, "N = %d (snapshot #%5d) max=%5.3f",
01426                 n_stars, init_flag + 1, lmax3d);
01427         lux_draw_image_string(win, 0.0, lmax3d*FAC3D, 0.5, temp_buffer, 0);
01428 
01429     } else {
01430 
01431         sprintf(temp_buffer, "N = %d  (snapshot #%5d)", n_stars, init_flag + 1);
01432         lux_draw_image_string(win, (xmin+xmax)/2.0, ymax, 0.5, temp_buffer, 0);
01433 
01434     }
01435 
01436     // Update the display.
01437 
01438     update_with_delay(win, delay_time);
01439 
01440     // Deal with user input.
01441 
01442     check_for_input(win, b, b_flag, f_flag);
01443 }
01444 
01445 /*-----------------------------------------------------------------------------
01446  *  main  --  driver to use  xstarplot()  as a tool. 
01447  *               The argument -a is interpreted as the axis along which to
01448  *            view the N-body system: the number of the coordinate axis along
01449  *            which the projection is directed, with {k,x,y} having a
01450  *            right-handend orientation.
01451  *               If an argument -l is provided, it defines the maximum
01452  *            lengths along the remaining axes.
01453  *-----------------------------------------------------------------------------
01454  */
01455 
01456 main(int argc, char** argv)
01457 {
01458     sdyn3 *b;
01459 
01460     // Establish some defaults:
01461 
01462     int   k = 3;                // Note x = 1, y = 2, z = 3 here!
01463     int   d = 2;
01464     float D = 0.0;
01465     int   cenergy = 0;
01466     float lmax = -1.0;
01467     int   point_mode = 0;
01468     float rel_point_size = -1.0;
01469     float scale = 1.0;
01470 
01471     bool  b_flag = TRUE;        // if TRUE, the background is black
01472     bool  f_flag = TRUE;        // if TRUE, the star is solid
01473     bool  o_flag = FALSE;       // if TRUE, output stdin to stdout
01474     bool  t_flag = FALSE;       // if TRUE, show tree structure
01475 
01476     extern char *poptarg;
01477     char* params = "a:d:D:efl:op:P:rs:t";
01478     int   c;
01479 
01480     while ((c = pgetopt(argc, argv, params)) != -1)
01481         switch(c) {
01482 
01483             case 'a': k = atoi(poptarg);        // projection axis [z]
01484                       break;
01485             case 'd': d = atoi(poptarg);        // number of dimensions [2]
01486                       break;
01487             case 'D': D = atof(poptarg);        // delay between frames [0]
01488                       break;
01489             case 'e': cenergy = 1;              // color by energy [no]
01490                       break;
01491             case 'f': f_flag = FALSE;           // fill star points [yes]
01492                       break;
01493             case 'l': lmax = atof(poptarg);     // axes +/- lmax [get from data]
01494                       break;
01495             case 'o': o_flag = TRUE;            // output stdin to stdout [no]
01496                       break;
01497             case 'p': point_mode = atoi(poptarg);       // point scale mode [0]
01498                       break;
01499             case 'P': rel_point_size = atof(poptarg);   // point scale factor
01500                       break;
01501             case 'r': b_flag = FALSE;           // black background [yes]
01502                       break;
01503             case 's': scale = atof(poptarg);    // rescale display [1.0]
01504                       break;
01505             case 't': t_flag = TRUE;            // show tree links [no]
01506                       break;
01507             case '?': params_to_usage(cerr, argv[0], params);
01508                       exit(0);
01509         }            
01510 
01511     // Reinterpret scaling if size reflects radius.
01512 
01513     if (point_mode == 2 && rel_point_size > 0)
01514         rel_point_size = 1.0/rel_point_size;
01515         
01516     int   gfx_counter = 0;
01517     while (b = get_sdyn3(cin)) {
01518         convert_relative_to_absolute(b);
01519         xstarplot(b, scale, k, d, lmax,
01520                   point_mode, rel_point_size, D, cenergy,
01521                   b_flag, f_flag, t_flag, gfx_counter++);
01522         if (o_flag) put_node(cout, *b);
01523         rmtree(b);
01524     }
01525 
01526     if (o_flag) cout << "End of data\n" << flush;
01527 
01528     // idle, wait to leave
01529 
01530     xstarplot(b, scale, k, d, lmax,
01531               point_mode, rel_point_size, D, cenergy,
01532               b_flag, f_flag, t_flag, -1);
01533 }