home *** CD-ROM | disk | FTP | other *** search
/ Ray Tracing Box / RAY_CD.mdf / RAYTRACE / RTAG / ORBIT.RTA < prev    next >
Encoding:
Text File  |  1993-03-03  |  3.3 KB  |  103 lines
  1. /*
  2.  *  RTAG source file for Orbital Simulation.
  3.  *  Created by Phillip H. Sherrod
  4.  */
  5.  
  6. //  ---  File declarations  ---
  7. bfile "orbit";            // Batch file is ORBIT.BAT
  8. //  ---  General simulation data  ---
  9. var numframes = 200;
  10. var endtime = 200;
  11. var timestep = .1;
  12. var sf = 5;            // Speed/distance factor
  13. //  ---  Planet orbital data  ---
  14. //  Planet 1
  15. var smaj1 = 6.443832;        // Semi-major axis
  16. var smin1 = 5.780208;        // Semi-minor axis
  17. var phase1 = 0;            // Phase angle
  18. //  Planet 2
  19. var smaj2 = 9.028421;
  20. var smin2 = 8.516810;
  21. var phase2 = 0;
  22. //  Planet 3
  23. var smaj3 = 12.960591;
  24. var smin3 = 11.520525;
  25. var phase3 = 0;
  26. //  Moon around planet 3
  27. var mr = 2;            // Radius of moon's orbit
  28. var phasem = 60;        // Orbital position angle when time = 0.
  29. var moonfreq = 5;        // Number of times moon orbits
  30. //  ---  Other variables  ---
  31. var x1,z1,f1,s1,d1;
  32. var x2,z2,f2,s2,d2;
  33. var x3,z3,f3,s3,d3;
  34. var xm,zm,sm;
  35. var time,frametime,framestep;
  36.  
  37. //  ---  Initialization  ---
  38. //  Compute position of foci (along x axis)
  39. f1 = sqrt(smaj1^2 - smin1^2);
  40. f2 = sqrt(smaj2^2 - smin2^2);
  41. f3 = sqrt(smaj3^2 - smin3^2);
  42. //  Angular velocity of moon
  43. sm = (360 * moonfreq) / endtime;
  44. //  Compute time between frames
  45. framestep = endtime / numframes;
  46. //  ---  Simulation procedure  ---
  47. for (time=0; time<endtime; time+=timestep) {
  48.     //  Compute current position
  49.     x1 = f1 + smaj1 * cos(phase1);
  50.     z1 = smin1 * sin(phase1);
  51.     x2 = f2 + smaj2 * cos(phase2);
  52.     z2 = smin2 * sin(phase2);
  53.     x3 = f3 + smaj3 * cos(phase3);
  54.     z3 = smin3 * sin(phase3);
  55.     xm = x3 + mr * cos(phasem);
  56.     zm = z3 + mr * sin(phasem);
  57.     //  See if it is time to output a frame
  58.     if (time >= frametime) {
  59.         nextframe;
  60.         //  Display our positions
  61.         printf("frame %3.0lf, time %5.1lf (%7.3lf,%7.3lf) (%7.3lf,%7.3lf) (%7.3lf,%7.3lf)\n",
  62.             curframe,time,x1,z1,x2,z2,x3,z3);
  63.         //  Generate batch file commands to declare position for frame.
  64.         bwrite("echo #declare x1 = `x1` > orbit.inc\n");
  65.         bwrite("echo #declare z1 = `z1` >> orbit.inc\n");
  66.         bwrite("echo #declare x2 = `x2` >> orbit.inc\n");
  67.         bwrite("echo #declare z2 = `z2` >> orbit.inc\n");
  68.         bwrite("echo #declare x3 = `x3` >> orbit.inc\n");
  69.         bwrite("echo #declare z3 = `z3` >> orbit.inc\n");
  70.         bwrite("echo #declare xm = `xm` >> orbit.inc\n");
  71.         bwrite("echo #declare zm = `zm` >> orbit.inc\n");
  72.         //  Generate batch file command to render this frame.
  73.         bwrite("call render orbit orbit`###`\n");
  74.         //  Remember when to generate next frame
  75.         frametime = frametime + framestep;
  76.     }
  77.     //  Compute distance from focus (the focus is at 0,0)
  78.     d1 = sqrt(x1*x1 + z1*z1);
  79.     d2 = sqrt(x2*x2 + z2*z2);
  80.     d3 = sqrt(x3*x3 + z3*z3);
  81.     //  Compute linear speed based on distance from focus (Kepler's law)
  82.     s1 = atan(1 / d1);
  83.     s2 = atan(1 / d2);
  84.     s3 = atan(1 / d3);
  85.     //  Since we are generating the ellipse from the center rather
  86.     //  than the focus, we must compute the angular speed necessary
  87.      //  to produce the desired linear speed.
  88.     d1 = sqrt((x1-f1)^2 + z1*z1);
  89.     d2 = sqrt((x2-f2)^2 + z2*z2);
  90.     d3 = sqrt((x3-f3)^2 + z3*z3);
  91.     s1 = (sf * s1) / d1;
  92.     s2 = (sf * s2) / d2;
  93.     s3 = (sf * s3) / d3;
  94.     //  Advance angles
  95.     phase1 = mod(phase1 + timestep * s1, 360);
  96.     phase2 = mod(phase2 + timestep * s2, 360);
  97.     phase3 = mod(phase3 + timestep * s3, 360);
  98.     phasem = mod(phasem + timestep * sm, 360);
  99. }
  100. //  Put last command in batch file.
  101. epilog;
  102. bwrite("call dodta ORBIT /S7\n");
  103.