attachment:tor08.c of Brown2013

   1 // tor08.c
   2 // v0.1 KHL 2009-04-21
   3 // gcc -o tor08 tor08.c -lm -lgd ; ./tor08
   4 // static drawing, though slow to compute
   5 // earth and toroid around orbit
   6 
   7 #define NAME     "tor08"
   8 #define NORB1     5
   9 #define NORB2     9
  10 #define DELTANG 0.2
  11 #define ORBWIDTH  5
  12 #define ANG0      0
  13 #define ANG1     90
  14 #define ANG2    142
  15 #define NELIP    32
  16 
  17 #include "tor00.hc"
  18 
  19 int main () {
  20    double  t_ang ;
  21    int     e_cnt ;
  22    int     orb ;
  23    int     norb ;
  24  
  25    gdPoint ellipse[NELIP+2];
  26 
  27    singlestart(NAME);
  28    framestart(    90, "M288, center + 4 toroid orbits");
  29    torstart();
  30 
  31    // draw back portion of the orbit toroid
  32   
  33    double acolor =  0.0 ;
  34    norb = NORB1 ;
  35    for( t_ang=ANG0 ; t_ang < ANG2 ; t_ang += DELTANG ) {
  36 
  37       double ecolor ;
  38       if( acolor < 0.5 ) { ecolor = gray ; }
  39       else {
  40          ecolor = lgray ;
  41          if( acolor > 14.5 ) { acolor -= 15.0 ; }
  42       }
  43       acolor += DELTANG ;
  44 
  45       double angle_rad=DEG2RAD( t_ang );
  46 
  47       // ellipse disk - since it is at an angle, we will have to draw it
  48       // as a filled polygon
  49     
  50       // ellipse disk on left
  51       for( e_cnt=0 ; e_cnt < NELIP ; e_cnt++ ) {
  52          double toa = angle_rad + (2.0*PI*e_cnt)/NELIP  ;
  53          double tos = td[ norb-1 ] * step ;                 // step size
  54          double tor = tos * cos( toa );
  55          double toy = tos * sin( toa );
  56          double tox = -(major+tor)*sin( angle_rad );
  57 	 double toz = -(major+tor)*cos( angle_rad );
  58 	 ellipse[e_cnt] = sr0( tox, toy, toz );
  59       }
  60       gdImagePolygon( im, ellipse, NELIP, ecolor );
  61 
  62       // ellipse disk on right
  63       for( e_cnt=0 ; e_cnt < NELIP ; e_cnt++ ) {
  64          double toa = -angle_rad + (2.0*PI*e_cnt)/NELIP  ;
  65          double tos = td[ norb-1 ] * step ;                 // step size
  66          double tor = tos * cos( toa );
  67          double toy = tos * sin( toa );
  68          double tox = -(major+tor)*sin( -angle_rad );
  69          double toz = -(major+tor)*cos( -angle_rad );
  70          ellipse[e_cnt] = sr0( tox, toy, toz );
  71       }
  72       gdImagePolygon( im, ellipse, NELIP, ecolor );
  73 
  74 
  75       for( orb=0 ; orb < norb ; orb++ ) {
  76          // orbit spot on left
  77          int    tcolor = tc[orb]; 
  78          double toa = angle_rad +  DEG2RAD( ta[orb] ) ;
  79          double tos = step * td[ orb ] ;                 // step size
  80          double tor = tos * cos( toa );
  81          double toy = tos * sin( toa );
  82          double tox = -(major+tor)*sin( angle_rad );
  83          double toz = -(major+tor)*cos( angle_rad );
  84          gdPoint pto = sr0( tox, toy, toz );
  85          gdImageFilledEllipse( im, pto.x, pto.y, ORBWIDTH, ORBWIDTH, tcolor);
  86          
  87          // orbit spot on right
  88          toa = -angle_rad +  DEG2RAD( ta[orb] ) ;
  89          tos = step * td[ orb ] ;                 // step size
  90          tor = tos * cos( toa );
  91          toy = tos * sin( toa );
  92          tox = -(major+tor)*sin( -angle_rad );
  93          toz = -(major+tor)*cos( -angle_rad );
  94          pto = sr0( tox, toy, toz );
  95          gdImageFilledEllipse( im, pto.x, pto.y, ORBWIDTH, ORBWIDTH, tcolor);
  96       }  
  97    }      
  98    drawearth();
  99 
 100    // draw the rest of the center orbit
 101    for( t_ang=ANG2 ; t_ang < 360.0-ANG2 ; t_ang += DELTANG ) {
 102       double  angle_rad=DEG2RAD( t_ang );
 103       double  tox = -major*sin( angle_rad );
 104       double  toz = -major*cos( angle_rad );
 105       gdPoint pto = sr0( tox, 0.0, toz );
 106       gdImageFilledEllipse( im, pto.x, pto.y, ORBWIDTH, ORBWIDTH, tc[0]);
 107    }
 108 
 109    frameend();
 110    return 0;
 111 }