head	1.15;
access;
symbols;
locks
	ids:1.15; strict;
comment	@ * @;


1.15
date	99.05.12.20.17.16;	author ids;	state Exp;
branches;
next	1.14;

1.14
date	99.04.27.18.50.08;	author ids;	state Exp;
branches;
next	1.13;

1.13
date	96.01.16.15.29.50;	author ids;	state Exp;
branches;
next	1.12;

1.12
date	93.08.05.20.45.18;	author ids;	state Exp;
branches;
next	1.11;

1.11
date	93.07.16.13.07.13;	author ids;	state Exp;
branches;
next	1.10;

1.10
date	93.07.15.21.59.57;	author ids;	state Exp;
branches;
next	1.9;

1.9
date	93.07.15.20.14.15;	author ids;	state Exp;
branches;
next	1.8;

1.8
date	93.07.15.17.51.09;	author ids;	state Exp;
branches;
next	1.7;

1.7
date	93.07.04.18.25.31;	author ids;	state Exp;
branches;
next	1.6;

1.6
date	93.06.26.20.57.05;	author ids;	state Exp;
branches;
next	1.5;

1.5
date	93.06.23.12.37.26;	author ids;	state Exp;
branches;
next	1.4;

1.4
date	93.06.22.13.29.32;	author ids;	state Exp;
branches;
next	1.3;

1.3
date	93.06.22.12.46.38;	author ids;	state Exp;
branches;
next	1.2;

1.2
date	93.06.22.08.35.58;	author ids;	state Exp;
branches;
next	1.1;

1.1
date	93.06.21.20.55.00;	author ids;	state Exp;
branches;
next	;


desc
@@


1.15
log
@at last! we seem to be able to dispense with the rigmarole of
 "plus give proper ANSI headers to..." as a hanger-on of "#include <math.h>".
 (It was things like sqrt, pow, sin, cos, atan2 that used to cause all the
  trouble - presumably they once had non-ANSI headers by default or something.)

Note: we do still have to define pi by hand. There's always something!
@
text
@#include "GraphicsForGravity.h"
#include "gravity_support.h"
#include "GraphicsSupport.h"
#include <stdio.h>
#include <ctype.h>
#include <math.h>
/* plus define pi */
#define pi		3.14159265358979323846
#include <X11/Xlib.h>	/* Xgraf unfortunately no longer does this for you! */
/*
 *	There are local copies of the "Xgraf and allied" files
 *	here in this directory
 *	(modified to use "filename.h" syntax
 *	rather than <filename.h> syntax as needed)
 *	because it would be far too complicated to try to
 *	resurrect and re-install
 *	the various mixed Modula-2 / C {compiling, linking, library-generating}
 *	mechanisms of old.
 *	Hence the "filename.h" syntax rather than <filename.h> syntax below.
 */
#include "Xgraf.h"
#include "Xcolours.h"

#define bool	int
#define TRUE	1
#define FALSE	0

#define WindowEventMask		\
		(MaskKeypress | MaskAllMotion | MaskMouseDown | MaskResize)
/* note: MaskAllMotion is for mouse movements */
#define WindowLeftX		(-1)
#define WindowTopY		(-1)
#define InitWindowWidth		500
#define InitWindowHeight	500
#define WindowDistance		500	/* eye-to-window, in pixels */

#define SMALL_FACTOR		0.1
/* for re-initializing the user translation distance and boost speed */

#define COARSE_ADJ_FACTOR	2.0
/* for changing the user translation distance and boost speed */

#define RotationAngle		(pi/128.0)


static Window		Win;
static int		WindowWidth, WindowHeight;
#define WindowShortSideLength	\
			(WindowWidth<WindowHeight ? WindowWidth : WindowHeight)
static GRAPHICS_STATE	gstate;



/* private forward references */
static void SetUpSensibleGraphicsState (STATE *s, GRAPHICS_STATE *gs);





void InitializeGraphics (STATE *AState)
{
  /* get window ready */
  InitialiseGraf();
  AllowNoEvent();	/* for clean semantics of EventPending() & GetEvent() */
  SetEventMask(WindowEventMask);
  SetWBackingStore();	/* asks server to try to provide off-screen pixmap */

  WindowWidth  = InitWindowWidth;
  WindowHeight = InitWindowHeight;

  Win = MakeWindow(RootWin(), "gravity", WindowLeftX, WindowTopY,
                   WindowWidth, WindowHeight);
  /*
  now set this window to be inverse video.
  NOTE: must do this here, after the call to MakeWindow; any earlier and the
        current window is still RootWin()!
  */
  SetColour(White);
  SetBackColour(Black);

  /* now get graphics state ready. */
  SetUpSensibleGraphicsState(AState, &gstate);
}





void DrawState (STATE *AState)
{
  /* MouseX and MouseY are for squirreling away the mouse position. */
  static int MouseX = 0.5*InitWindowWidth;
  static int MouseY = 0.5*InitWindowHeight;

  /*
  UserDimension is for storing a user dimension number from 1 to NDIMS, and is
  initialized to 0 to show it hasn't yet been obtained from the keyboard.
  */
  static int UserDimension = 0;

  /* TrailsWanted stores whether or not the bodies are to leave trails. */
  static bool TrailsWanted = TRUE;

  Window     dummywin;
  EventType  event;
  EventData  eventdata;

  GRAPHICS_VECTOR apparent_user_velocity;	/* in user coordinates */

  GRAPHICS_VECTOR boost_velocity;		/* in actual coordinates */
  GRAPHICS_VECTOR apparent_boost_velocity;	/* in user coordinates */

  VECTOR P;		/* a body in user coordinates */
  double     Rx, Ry;	/* 2-D retinal coordinates */
  double     Rlength;	/* used for moving towards a given retinal direction */
  int        i,j,k,d;
  GRAPHICS_MATRIX temp_matrix;


  /*
  first update the graphics state to take into account the drift of
  user position since the last time we updated it.
  */
  MultiplyGraphicsVectorByMatrix
   (gstate.user_velocity, gstate.M, apparent_user_velocity);
  for (d=0; d<NDIMS+1; d++) {
    gstate.M[ NDIMS ][ d ] -=
     apparent_user_velocity[d] * (AState->Time - gstate.TimeLastUpdated);
  }/*endfor*/
  gstate.TimeLastUpdated = AState->Time;


  /*
  now see if any events need special attention. Note that EventExpose,
  perhaps surprisingly, does not need special attention: we're going to
  draw the latest state anyway.
  (albeit without clearing the window, if trails are wanted)
  */
  while (EventPending()) {
    GetEvent(&dummywin, &event, &eventdata);
    switch (event) {
      case EventMove :	/* mouse movement, not moving of the whole window! */

        /* squirrel away the mouse position for future use by "m" or "M" */
        MouseX = eventdata.x;
        MouseY = eventdata.y;
        break;

      case EventResize :

        WindowWidth  = eventdata.w;
        WindowHeight = eventdata.h;
        ClearWindow();
        break;

      case EventKeypress :

        #define FLIP_FACTOR	( islower(eventdata.ch) ? 1.0 : -1.0 )

        switch (eventdata.ch) {
          case '0' : case '1' : case '2' : case '3' : case '4' :
          case '5' : case '6' : case '7' : case '8' : case '9' :
            UserDimension = eventdata.ch == '0' ? 10 : eventdata.ch - '0';
            if (UserDimension > NDIMS) {
              beep();
              UserDimension = 0;
            }/*endif*/
            break;

          case 'b' : case 'B' :
            if (UserDimension > 0) {
              for (d=0; d<NDIMS+1; d++) {
                apparent_boost_velocity[d] =
                 d==UserDimension-1 ? FLIP_FACTOR * gstate.BoostSpeed : 0.0;
              }/*endfor*/
              InvertOrthonormal(gstate.M, temp_matrix);
              MultiplyGraphicsVectorByMatrix
               (apparent_boost_velocity, temp_matrix, boost_velocity);
              AddGraphicsVectors
               (gstate.user_velocity, boost_velocity, gstate.user_velocity);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 'c'            : ClearWindow();
                                break;

          case 'C'            : TrailsWanted = ! TrailsWanted;
                                break;

          case 'g' : case 'G' : /* go to a good vantage point */
                                SetUpSensibleGraphicsState(AState, &gstate);
                                ClearWindow();
                                break;

          case 'h' : case 'H' :
            printf("*****HELP MENU*****\n");
            printf("h,H	display this help menu\n");
            printf("q,Q	quit\n\n");
            printf("In the examples below any valid dimension-numbers can be used.\n");
            if (NDIMS == 10) {
              printf("[For dimension 10 press 0.]\n");
            }/*endif*/
            printf("Also, t,T,r,R,b,B will auto-repeat without repeating the dimension number.\n\n");
            printf("m	translate in direction indicated by mouse cursor\n");
            printf("M	translate in opposite direction to mouse cursor\n");
            printf("3t	translate along positive 3-axis\n");
            printf("3T	translate along negative 3-axis\n");
            printf("LeftMouseButton: rotate to face direction clicked on\n");
            printf("2r	rotate towards positive 2-axis\n");
            printf("2R	rotate towards negative 2-axis\n");
            printf("3b	boost along positive 3-axis\n");
            printf("3B	boost along negative 3-axis\n");
            printf("c  	clear display window\n");
            printf("C  	toggle between clearing vs. leaving trails\n");
            printf("x	reduce translation distance\n");
            printf("X	increase translation distance\n");
            printf("v	reduce boost speed\n");
            printf("V	increase boost speed\n");
            printf("g,G	go to a good vantage point\n\n\n");
            break;

          case 'm' : case 'M' :
            Rx =   (MouseX + 0.5 - 0.5*WindowWidth ) / WindowDistance;
            Ry = - (MouseY + 0.5 - 0.5*WindowHeight) / WindowDistance;
            Rlength = sqrt(Rx*Rx+Ry*Ry+1.0);

            /* now move the correct distance in the required direction */
            gstate.M[ NDIMS ][0] -=
             FLIP_FACTOR *  Rx / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][1] -=
             FLIP_FACTOR *  Ry / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][ NDIMS-1 ] -=
             FLIP_FACTOR * 1.0 / Rlength * gstate.TranslationDistance;

            ClearWindow();
            break;

          case 'q' : case 'Q' : CloseDownGraf();
                                exit(0);
                                break;

          case 'r' : case 'R' :
            if (UserDimension > 0 && UserDimension < NDIMS) {
              CreateRotationMatrix
               (NDIMS-1, UserDimension-1, - FLIP_FACTOR * RotationAngle,
                temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 't' : case 'T' :
            if (UserDimension > 0) {
              CreateTranslationMatrix
               (UserDimension-1, - FLIP_FACTOR * gstate.TranslationDistance,
                temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 'x'            : gstate.TranslationDistance /= COARSE_ADJ_FACTOR;
                                break;

          case 'X'            : gstate.TranslationDistance *= COARSE_ADJ_FACTOR;
                                break;

          case 'v'            : gstate.BoostSpeed /= COARSE_ADJ_FACTOR;
                                break;

          case 'V'            : gstate.BoostSpeed *= COARSE_ADJ_FACTOR;
                                break;

          default             : break;
        }/*endswitch*/
        break;

      case EventMouseDown :

        switch (eventdata.button) {
          case 1 :
            Rx =   (eventdata.x + 0.5 - 0.5*WindowWidth ) / WindowDistance;
            Ry = - (eventdata.y + 0.5 - 0.5*WindowHeight) / WindowDistance;
            if (Rx != 0.0 || Ry != 0.0) {
              /* first twist so retinal point is on positive 1-axis */
              CreateRotationMatrix
               (1, 0, atan2(Ry,Rx), temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);

              /* now rotate point to exact centre of view */
              CreateRotationMatrix
               (0, NDIMS-1, atan2(sqrt(Rx*Rx+Ry*Ry), 1.0), temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);

              /* now undo twist above. */
              CreateRotationMatrix
               (0, 1, atan2(Ry,Rx), temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);

              ClearWindow();
            }/*endif*/
            break;

          default :
            break;
        }/*endswitch*/
        break;

      default :

        break;

    }/*endswitch*/
  }/*endwhile*/


  if (! TrailsWanted) {
    ClearWindow();
  }/*endif*/


  /* now, just draw all the bodies! */
  for (i=0; i<AState->n; i++) {
    /*
    turn body i into user coordinates. We use a tight little loop here as
    it's actually clearer to do this than to turn everything into graphics
    types and then use the graphics matrix/vector functions.
    (not to mention faster, in this "fairly inner" position.....)
    */
    CreateZeroVector(P);
    for (j=0; j<NDIMS; j++) {
      for (k=0; k<NDIMS+1; k++) {
        P[j] +=
         (k<NDIMS ? AState->bodydata[i].position[k] : 1.0) * gstate.M[k][j];
      }/*endfor*/
    }/*endfor*/

    /* check if point is in front of user */
    if (P[NDIMS-1] > 0) {
      /* find retinal coordinates */
      Rx = P[0] / P[NDIMS-1];
      Ry = P[1] / P[NDIMS-1];

      /* now plot in Xgraf coordinates */
      DrawPoint( (int) (0.5*WindowWidth  + Rx*WindowDistance),
                 (int) (0.5*WindowHeight - Ry*WindowDistance)
               );
    }/*endif*/
  }/*endfor*/
}





static void SetUpSensibleGraphicsState
             (STATE *AState, GRAPHICS_STATE *AGraphicsState)
/*
Treats its first arg as read-only.
Sets up sensible values in the indicated graphics state, with the user
co-moving with the system of bodies at a well-chosen vantage point.
*/
{
  double total_mass;
  VECTOR total_moment;
  VECTOR total_momentum;
  VECTOR centre_of_mass_position;
  VECTOR centre_of_mass_velocity;
  VECTOR position_from_CM, velocity_from_CM;
  double max_CM_disp, max_CM_speed;
  int    i, d;

  /*
  find centre-of-mass position and velocity, and maximum deviations therefrom,
  in order to set up user to be co-moving with the system
  at a well-chosen vantage point.
  */
  total_mass   = 0.0;
  CreateZeroVector(total_moment);
  CreateZeroVector(total_momentum);
  max_CM_disp  = 0.0;
  max_CM_speed = 0.0;

  for (i=0; i<AState->n; i++) {
    total_mass += AState->bodydata[i].mass;
    for (d=0; d<NDIMS; d++) {
      total_moment[d]   += AState->bodydata[i].mass * AState->bodydata[i].position[d];
      total_momentum[d] += AState->bodydata[i].mass * AState->bodydata[i].velocity[d];
    }/*endfor*/;
  }/*endfor*/

  for (d=0; d<NDIMS; d++) {
    centre_of_mass_position[d] = total_moment[d] / total_mass;
    centre_of_mass_velocity[d] = total_momentum[d] / total_mass;
  }/*endfor*/;

  for (i=0; i<AState->n; i++) {
    for (d=0; d<NDIMS; d++) {
      position_from_CM[d] =
       AState->bodydata[i].position[d] - centre_of_mass_position[d];
      velocity_from_CM[d] =
       AState->bodydata[i].velocity[d] - centre_of_mass_velocity[d];
    }/*endfor*/;

    if (sqrt(lengthsquared(position_from_CM)) > max_CM_disp) {
      max_CM_disp = sqrt(lengthsquared(position_from_CM));
    }/*endif*/

    if (sqrt(lengthsquared(velocity_from_CM)) > max_CM_speed) {
      max_CM_speed = sqrt(lengthsquared(velocity_from_CM));
    }/*endif*/
  }/*endfor*/;


  /* now use this hard-won data to set a sensible initial user position etc. */

  /* first put user exactly *on* the centre-of-mass..... */
  CreateIdentityMatrix(AGraphicsState->M);
  for(d=0; d<NDIMS+1; d++) {
    AGraphicsState->M[ NDIMS ][ d ] =
     d<NDIMS ? (- centre_of_mass_position[d]) : 1.0;
  }/*endfor*/

  /* .....and now back away a well-chosen length. */
  AGraphicsState->M[ NDIMS ][ NDIMS-1 ] +=
   sqrt(16.0 * pow(WindowDistance, 2.0) / pow(WindowShortSideLength, 2.0) + 4.0)
    * max_CM_disp;

  /* start user off co-moving with the system. */
  for (d=0; d<NDIMS+1; d++) {
    AGraphicsState->user_velocity[d] =
     d<NDIMS ? centre_of_mass_velocity[d] : 0.0;
  }/*endfor*/

  /* set up sensible user translation distance and boost speed */
  AGraphicsState->TranslationDistance = SMALL_FACTOR * max_CM_disp;
  AGraphicsState->BoostSpeed          = SMALL_FACTOR * max_CM_speed;

  /*
  finally, record the time, to enable future updates of the graphics state.
  */
  AGraphicsState->TimeLastUpdated = AState->Time;
}
@


1.14
log
@hacked to cope with the local copies of the "Xgraf and allied" files
 which we unfortunately need here in this directory for now:

> /*
>  *    There are local copies of the "Xgraf and allied" files
>  *    here in this directory
>  *    (modified to use "filename.h" syntax
>  *    rather than <filename.h> syntax as needed)
>  *    because it would be far too complicated to try to
>  *    resurrect and re-install
>  *    the various mixed Modula-2 / C {compiling, linking, library-generating}
>  *    mechanisms of old.
>  *    Hence the "filename.h" syntax rather than <filename.h> syntax below.
>  */

If and when the Xgraf and allied files get properly re-installed
 in something like the old style, we can presumably change the syntax back.
@
text
@d7 1
a7 4
/* plus give proper ANSI header to sqrt, pow, atan2, and define pi */
extern double sqrt (double);
extern double pow  (double, double);
extern double atan2(double, double);
@


1.13
log
@added explicit include of X11 stuff (Xgraf unfortunately no longer does
 this for you!).
@
text
@d13 13
a25 2
#include <Xgraf.h>
#include <Xcolours.h>
@


1.12
log
@added CloseDownGraf() : it does matter! (to release server memory, etc.)
@
text
@d12 1
@


1.11
log
@cleaner, shorter code by consolidating code for inverse pairs of events
 ('m', 'M', 't', 'T', etc.) into one piece of code for each pair.
better order in help menu
C tidyups.
@
text
@d233 2
a234 1
          case 'q' : case 'Q' : exit(0);
@


1.10
log
@added facility for moving in the direction indicated with the mouse cursor.
@
text
@d4 2
d150 2
d162 1
a162 18
          case 'b' :
            if (UserDimension > 0) {
              for (d=0; d<NDIMS+1; d++) {
                apparent_boost_velocity[d] =
                 d==UserDimension-1 ? gstate.BoostSpeed : 0.0;
              }/*endfor*/
              InvertOrthonormal(gstate.M, temp_matrix);
              MultiplyGraphicsVectorByMatrix
               (apparent_boost_velocity, temp_matrix, boost_velocity);
              AddGraphicsVectors
               (gstate.user_velocity, boost_velocity, gstate.user_velocity);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 'B' :
d166 1
a166 1
                 d==UserDimension-1 ? - gstate.BoostSpeed : 0.0;
d199 2
d203 1
a203 2
            printf("m	translate in direction indicated by mouse cursor\n");
            printf("M	translate in opposite direction to mouse cursor\n");
a205 1
            printf("LeftMouseButton: rotate to face direction clicked on\n");
d217 1
a217 14
          case 'm' :
            Rx =   (MouseX + 0.5 - 0.5*WindowWidth ) / WindowDistance;
            Ry = - (MouseY + 0.5 - 0.5*WindowHeight) / WindowDistance;
            Rlength = sqrt(Rx*Rx+Ry*Ry+1.0);

            /* now move the correct distance in the required direction */
            gstate.M[ NDIMS ][0] -= Rx / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][1] -= Ry / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][ NDIMS-1 ] -= 1.0 / Rlength * gstate.TranslationDistance;

            ClearWindow();
            break;

          case 'M' :
d223 6
a228 3
            gstate.M[ NDIMS ][0] += Rx / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][1] += Ry / Rlength * gstate.TranslationDistance;
            gstate.M[ NDIMS ][ NDIMS-1 ] += 1.0 / Rlength * gstate.TranslationDistance;
d236 1
a236 12
          case 'r' :
            if (UserDimension > 0 && UserDimension < NDIMS) {
              CreateRotationMatrix
               (NDIMS-1, UserDimension-1, - RotationAngle, temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 'R' :
d239 2
a240 12
               (NDIMS-1, UserDimension-1,   RotationAngle, temp_matrix);
              MultiplyMatrices(gstate.M, temp_matrix, gstate.M);
              ClearWindow();
            }else{
              beep();
            }/*endif*/
            break;

          case 't' :
            if (UserDimension > 0) {
              CreateTranslationMatrix
               (UserDimension-1, - gstate.TranslationDistance, temp_matrix);
d248 1
a248 1
          case 'T' :
d251 2
a252 1
               (UserDimension-1,   gstate.TranslationDistance, temp_matrix);
@


1.9
log
@added facility to turn immediately to face the direction clicked on with the
 left mouse button. (You can then travel straight ahead and visit whatever
 object you turned to face, if any.)
@
text
@d17 3
a19 1
#define WindowEventMask		(MaskKeypress | MaskMouseDown | MaskResize)
d81 4
d105 1
d132 7
d214 2
d218 1
a218 1
            printf("LeftMouseButton: rotate to face position clicked on\n");
d228 26
@


1.8
log
@added toggle facility to choose whether or not bodies are to leave trails.
@
text
@d5 4
a8 3
/* plus give proper ANSI header to sqrt and pow, and define pi */
extern double sqrt(double);
extern double pow (double, double);
d17 1
a17 1
#define WindowEventMask		(MaskKeypress | MaskResize)
d202 1
d274 31
@


1.7
log
@added resize capability, with sensible handling of viewpoint.
made inverse video (it's bodies in space after all, innit?!).
@
text
@d12 4
d84 3
d118 2
a119 1
  draw the latest state (albeit without clearing the window) anyway.
d177 4
a180 1
          case 'c' : case 'C' : ClearWindow();
d203 2
a204 1
            printf("c,C	clear display window\n");
d281 5
@


1.6
log
@new cleaner C conventions used, bearing a much stronger resemblance to
 Modula's modularization.
@
text
@d5 1
a5 1
/* plus give proper ANSI header to sqrt, and define pi */
d7 1
d10 1
d12 6
a17 4
#define WindowLeftX	(-1)
#define WindowTopY	(-1)
#define WindowSize	500
#define RetinalWidth	1.0
a18 3
#define BACK_AWAY_FACTOR	5.0
/* for backing away from centre-of-mass initially, to get a good view */

d28 6
a33 2
static Window         Win;
static GRAPHICS_STATE gstate;
d48 1
d50 4
d55 8
a62 1
                   WindowSize, WindowSize);
d64 1
a64 1
  /* get graphics state ready */
d116 7
d294 2
a295 2
      DrawPoint( (int) (( Rx/RetinalWidth+0.5) * WindowSize),
                 (int) ((-Ry/RetinalWidth+0.5) * WindowSize)
d374 3
a376 1
  AGraphicsState->M[ NDIMS ][ NDIMS-1 ] += BACK_AWAY_FACTOR * max_CM_disp;
@


1.5
log
@QUITE MAJOR REVISION:
 * added user control of graphics state, covering translating, rotating,
   boosting velocity, and controlling the size of the translates/boosts.
 * added help menu.
 * cleaned up computation of good graphics state, to compute max. velocities
   as well as positions, and to compute both relative to centre-of-mass,
   rather than the old span-of-system method.
@
text
@d31 2
d36 2
d281 7
a287 1
void SetUpSensibleGraphicsState (STATE *AState, GRAPHICS_STATE *AGraphicsState)
@


1.4
log
@added call to SetWBackingStore().
[Note that the server need not pay any attention to this request!]
@
text
@d3 1
d5 1
a5 1
/* plus give proper ANSI header to sqrt */
d7 1
a7 1
#define HASPROTOS
d15 1
a15 1
#define BACK_AWAY_FACTOR	3.0
d18 2
d21 6
d53 9
a61 7
  Window    dummywin;
  EventType event;
  EventData eventdata;

  VECTOR    P;		/* a body in user coordinates */
  double    Rx, Ry;	/* 2-D retinal coordinates */
  int       i,j,k,d;
d63 11
d78 2
d82 1
a82 1
     gstate.user_velocity[d] * (AState->Time - gstate.TimeLastUpdated);
d86 1
d98 47
a144 1
          case 'r' : case 'R' : /* recompute a good graphics state */
d149 77
a225 1
          case 'c' : case 'C' : ClearWindow();
d228 1
a228 1
          case 'q' : case 'Q' : exit(0);
d242 1
d245 6
a250 1
    /* turn body i into user coordinates */
d279 1
a279 1
  double total_mass, span_of_system;
d284 3
a286 2
  int    i, j, d;
  VECTOR i_to_j;
d289 1
a289 1
  find centre-of-mass position and velocity, and span of system,
d293 1
a293 2
  total_mass = 0.0;
  span_of_system = 0.0;
d296 2
a304 8
    for (j=i+1; j<AState->n; j++) {
      for (d=0; d<NDIMS; d++) {
        i_to_j[d] = AState->bodydata[j].position[d] - AState->bodydata[i].position[d];
      }/*endfor*/;
      if (sqrt(lengthsquared(i_to_j)) > span_of_system) {
        span_of_system = sqrt(lengthsquared(i_to_j));
      }/*endif*/
    }/*endfor*/;
d312 17
d333 4
a336 9
  for (i=0; i<NDIMS; i++) {
    for (j=0; j<NDIMS+1; j++) {
      AGraphicsState->M[i][j] = (i==j);
    }/*endfor*/
  }/*endfor*/

  for(j=0; j<NDIMS+1; j++) {
    AGraphicsState->M[ NDIMS ][ j ] =
     j<NDIMS ? (- centre_of_mass_position[j]) : 1.0;
d340 1
a340 1
  AGraphicsState->M[ NDIMS ][ NDIMS-1 ] += BACK_AWAY_FACTOR * span_of_system;
d347 4
@


1.3
log
@used type GRAPHICS_STATE
added keypress 'R' to recompute good vantage-point at any time
cleaner "proceduralization"
better comments.
@
text
@d30 1
@


1.2
log
@better procedure names, following the "IC DoC conventions"
(Put___ for standard output, Draw___ for graphics).
@
text
@d6 1
d18 2
a19 2
typedef double          GRAPHICS_VECTOR [NDIMS+1];
typedef GRAPHICS_VECTOR GRAPHICS_MATRIX [NDIMS+1];
a21 4
static Window          Win;
static GRAPHICS_MATRIX M;	/* actual-to-user, for time TimeLastUpdated */
static VECTOR          user_velocity;
static double          TimeLastUpdated;
a24 2


a26 8
  double total_mass, span_of_system;
  VECTOR total_moment;
  VECTOR total_momentum;
  VECTOR centre_of_mass_position;
  VECTOR centre_of_mass_velocity;
  int    i, j, d;
  VECTOR i_to_j;

d33 2
a34 55
  /*
  find centre-of-mass position and velocity, and span of system,
  in order to initialize user as drifting along with the system
  (same velocity), at a well-chosen vantage point.
  */
  total_mass = 0.0;
  span_of_system = 0.0;
  CreateZeroVector(total_moment);
  CreateZeroVector(total_momentum);

  for (i=0; i<AState->n; i++) {
    total_mass += AState->bodydata[i].mass;
    for (d=0; d<NDIMS; d++) {
      total_moment[d]   += AState->bodydata[i].mass * AState->bodydata[i].position[d];
      total_momentum[d] += AState->bodydata[i].mass * AState->bodydata[i].velocity[d];
    }/*endfor*/;
    for (j=i+1; j<AState->n; j++) {
      for (d=0; d<NDIMS; d++) {
        i_to_j[d] = AState->bodydata[j].position[d] - AState->bodydata[i].position[d];
      }/*endfor*/;
      if (sqrt(lengthsquared(i_to_j)) > span_of_system) {
        span_of_system = sqrt(lengthsquared(i_to_j));
      }/*endif*/
    }/*endfor*/;
  }/*endfor*/

  for (d=0; d<NDIMS; d++) {
    centre_of_mass_position[d] = total_moment[d] / total_mass;
    centre_of_mass_velocity[d] = total_momentum[d] / total_mass;
  }/*endfor*/;


  /* now use this hard-won data to set a sensible initial user position etc. */

  /* first put user exactly *on* the centre-of-mass..... */
  for (i=0; i<NDIMS; i++) {
    for (j=0; j<NDIMS+1; j++) {
      M[i][j] = (i==j);
    }/*endfor*/
  }/*endfor*/

  for(j=0; j<NDIMS+1; j++) {
    M[ NDIMS ][ j ] = j<NDIMS ? (- centre_of_mass_position[j]) : 1.0;
  }/*endfor*/

  /* .....and now back away a well-chosen length. */
  M[ NDIMS ][ NDIMS-1 ] += BACK_AWAY_FACTOR * span_of_system;

  /* start user off drifting at equal velocity with the system. */
  for (d=0; d<NDIMS; d++) {
    user_velocity[d] = centre_of_mass_velocity[d];
  }/*endfor*/

  /* finally, record the time, to enable future updates of the matrix M. */
  TimeLastUpdated = AState->Time;
d52 2
a53 2
  first update the matrix to take into account the drift of user position
  since the last time we updated it.
d55 3
a57 2
  for (d=0; d<NDIMS; d++) {
    M[ NDIMS ][ d ] -= user_velocity[d] * (AState->Time - TimeLastUpdated);
d59 1
a59 1
  TimeLastUpdated = AState->Time;
d72 5
d79 1
d82 1
d94 1
a94 1
  /* now, just draw the latest state! */
d100 2
a101 1
        P[j] += (k<NDIMS ? AState->bodydata[i].position[k] : 1.0) * M[k][j];
d117 75
@


1.1
log
@Initial revision
@
text
@d107 1
a107 1
void PlotGraphically (STATE *AState)
@