// ss6a.c
// compile:  gcc -o ss6a ss6a.c -lgd -lpng -lm
//
//  The perspective sucks
//  The light pressure of tilted thrusters is treated here as constant. 
//  The actual light pressure is proportional to cos^2( sun angle )
//  This assumes a 4:3 ratio, with gauges down the right side

#define  NAME  "ss6a"
#define  RATE   10  

#define FNT    "DejaVuMonoSans"
#define SCALE   0.7              // size of serversat relative to YSIZE

#define FSZ     26               // this scales the font and the whole drawing

#define Z0      24.0             // depth scaling for cheesy perspective
#define ACC      8.0             // acceleration in degrees/minute^2
#define TSTART -60.0             //
#define TIME1  150.0             // acceleration time in seconds
#define TIME2  (2.0*TIME1)       // deceleration time
#define TEND   360.0             //
#define TDELTA   2.0             // display timestep
#define NSTARS 500               // background stars (this shows transparency)
#define MOVE   0.3               // how much the background stars move
#define THICK    5               // dummy thickness in pixels
#define FILTER 0.3               //

#define LINETHICK 1
#define XSIZE   1024
#define YSIZE   768
#define RIGHT   (YSIZE-20)
#define MINUTE  60.0
#define HOUR    3600.0
#define MINUTES (TIME1/MINUTE)

#include "gd.h"
#include <math.h>
#include <string.h>
#include <stdlib.h>
#include <sys/stat.h>
#include <sys/types.h>


#define PI 3.14159265359
#define DEG2RAD(x) ((x)*PI/180.)

// #define WIRING 1              // show solar cell wires

// ----------------------------------------------------
char     rmmkdir[200] ;
char     pngfmt[ 200] ;
char     png2swf[200] ;

gdImagePtr im;
int      black, white, gray   ;
int      lred, lgreen         ;
int      dgray, trans         ;
int      xgauge=(XSIZE+YSIZE)/2 ;// down the right edge
double   scale                ; // scale serversat in pixels

int      spdcolor             ;
int      yellow[256]          ; // Yellow color range
// starting values
int      rate      = RATE     ;

double   ang  = 0.0           ; //  angle in radians
double   angc                 ; //  cos(ang)
double   angs                 ; //  sin(ang)
double   angt                 ; //  tan(ang)
 
//          | bottom |  front |   top  |  back  |
// an2    -135      -45      +45     +135     -135
// an2a          -        +        +        -
// an2b          +        +        -        -

double   an2                  ; //  double ang
double   pi4                  ; //  
double   pi                   ; //  
double   an2c                 ; //  cos(an2)
double   an2s                 ; //  sin(an2)
double   an2a                 ; //  an2 test a
double   an2b                 ; //  an2 test b
 
double   time      = 0.0      ;
double   angle_vel = 0.0      ;
double   angle_acc = 0.0      ;
double   angle     = 0.0      ;

double   min_vel   = 0.0      ;
double   max_vel   = ACC*TIME1/MINUTE ;
double   angle_min = 0.0      ;
double   angle_max = ACC*TIME1*TIME1/(MINUTE*MINUTE);
double   min_acc   = -ACC     ;
double   max_acc   = ACC      ;

double   thick     = THICK ;

//--------------------------------------------------------------------
// fixed rectangle subroutines, swap order if needed 
//  gdImageRectangle silently fails unless it gets left, bottom, right, top

void RectE( gdImagePtr imm,
            double rx1, double ry1, double rx2, double ry2, int rcolor ) {
   if( rx1 < rx2 ) {
      if( ry1 < ry2 ) {
         gdImageRectangle( imm, rx1, ry1, rx2, ry2, rcolor );
      } else {
         gdImageRectangle( imm, rx1, ry2, rx2, ry1, rcolor );
      }
   } else {
      if( ry1 < ry2 ) {
         gdImageRectangle( imm, rx2, ry1, rx1, ry2, rcolor );
      } else {
         gdImageRectangle( imm, rx2, ry2, rx1, ry1, rcolor );
      }
   }
}

   
void RectF( gdImagePtr imm,
            double rx1, double ry1, double rx2, double ry2, int rcolor ) {
   if( rx1 < rx2 ) {
      if( ry1 < ry2 ) {
         gdImageFilledRectangle( imm, rx1, ry1, rx2, ry2, rcolor );
      } else {
         gdImageFilledRectangle( imm, rx1, ry2, rx2, ry1, rcolor );
      }
   } else {
      if( ry1 < ry2 ) {
         gdImageFilledRectangle( imm, rx2, ry1, rx1, ry2, rcolor );
      } else {
         gdImageFilledRectangle( imm, rx2, ry2, rx1, ry1, rcolor );
      }
   }
}
   
// draw gauges ---------------------------------
// uses globals

#define  YCLOCK   (4*FSZ)
#define  YACC     (10*FSZ)
#define  YVEL     (16*FSZ)
#define  YANG     (22*FSZ)
#define  F2       (4*FSZ/5)
#define  F4       (FSZ/4)

#define  RCLOCK    (3*FSZ+8)
#define  SECHAND   (3*FSZ)
#define  MINHAND1  (2*FSZ)
#define  MINHAND2  (2*FSZ+5)
#define  MINWID    40

int      gleft  = YSIZE + 2*FSZ ;
int      gright = XSIZE - 2*FSZ ;


int  gscale( double gvar, double gmin, double gmax ) {
    return ( (int) ( gleft + (gvar-gmin)
                           * ( (double)(gright-gleft))/(gmax-gmin) ) ) ;
}

void drawgauges()  {

   double wangle ;   // working angle
   int    cx1, cy1, cx2, cy2 ;
   char   s[80];
   int    ic ;
   double ltime = time ;
   if( ltime < 0.0   ) { ltime = 0.0   ; }
   if( ltime > TIME2 ) { ltime = TIME2 ; }

   // clock at xgauge, 100 -----------------------------

   gdImageArc(  im, xgauge, YCLOCK, 2*RCLOCK, 2*RCLOCK, 0, 360, white );

   // clock / 60
   for( ic = 0;  ic < 60 ; ic += 1 ) {
      wangle  = PI * ic / 30 ;
      if ( ic%15 == 0 ) {
         cx1     = xgauge + (int)(0.5 + (RCLOCK-6)*sin(wangle)) ;
         cy1     = YCLOCK - (int)(0.5 + (RCLOCK-6)*cos(wangle)) ;
      }
      else if( ic%5 == 0 ) {
         cx1     = xgauge + (int)(0.5 + (RCLOCK-4)*sin(wangle)) ;
         cy1     = YCLOCK - (int)(0.5 + (RCLOCK-4)*cos(wangle)) ;
      }
      else {
         cx1     = xgauge + (int)(0.5 + (RCLOCK-2)*sin(wangle)) ;
         cy1     = YCLOCK - (int)(0.5 + (RCLOCK-2)*cos(wangle)) ;
      }
      cx2     = xgauge +  (int)(0.5 + RCLOCK   *sin(wangle) ) ;
      cy2     = YCLOCK -  (int)(0.5 + RCLOCK   *cos(wangle) ) ;
      gdImageLine( im, cx1, cy1, cx2, cy2, white ) ;
   }

   // second hand
   wangle  = 2.0 * PI * ltime / MINUTE ;
   cx1     = xgauge+(int)(0.5 + SECHAND*sin(wangle) ) ;
   cy1     = YCLOCK-(int)(0.5 + SECHAND*cos(wangle) ) ;
   gdImageLine( im, xgauge, YCLOCK, cx1, cy1, white ) ;


   // minute hand (two strokes )
   wangle  = 2.0 * PI * (ltime + MINWID) / HOUR ;
   cx1     = xgauge+(int)(0.5 + MINHAND1*sin(wangle) ) ;
   cy1     = YCLOCK-(int)(0.5 + MINHAND1*cos(wangle) ) ;

   wangle  = 2.0 * PI * (ltime         ) / HOUR ;
   cx2     = xgauge+(int)(0.5 + MINHAND2*sin(wangle) ) ;
   cy2     = YCLOCK-(int)(0.5 + MINHAND2*cos(wangle) ) ;
   gdImageLine( im, xgauge, YCLOCK, cx1, cy1, white ) ;
   gdImageLine( im,    cx2,    cy2, cx1, cy1, white );
   
   wangle  = 2.0 * PI * (ltime - MINWID) / HOUR ;
   cx1     = xgauge+(int)(0.5 + MINHAND1*sin(wangle) ) ;
   cy1     = YCLOCK-(int)(0.5 + MINHAND1*cos(wangle) ) ;
   gdImageLine( im,    cx2,    cy2, cx1, cy1, white );
   gdImageLine( im, xgauge, YCLOCK, cx1, cy1, white ) ;

   // text clock

   sprintf( s, "%3.0f secs", ltime );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, xgauge-2*FSZ, 9*FSZ, s );

   // angular acceleration gauge ----------------

   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0,
                              YSIZE, YACC+2*FSZ-F4, "acceleration" );
   gdImageStringFT( im, NULL, white, FNT, F2, 0.0,
                              YSIZE, YACC+5*FSZ, "degrees/min^2" );


   sprintf( s, "%2.0f", min_acc );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, YSIZE,  YACC+4*FSZ, s );

   sprintf( s, "%2.0f", max_acc );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, gright+F4, YACC+4*FSZ, s );
   
   cx1 =  gscale( angle_acc, min_acc, max_acc ) ;
   cx2 =  gscale(         0, min_acc, max_acc ) ;
   RectE( im, gleft, YACC+2*FSZ, gright, YACC+4*FSZ, white );
   RectF( im,   cx1, YACC+2*FSZ,    cx2, YACC+4*FSZ, white );

   // angular velocity gauge ---------------------

   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0,
                              YSIZE, YVEL+2*FSZ-F4, "velocity" );
   gdImageStringFT( im, NULL, white, FNT, F2, 0.0,
                              YSIZE, YVEL+5*FSZ, "degrees/min" );

   sprintf( s, "%2.0f", min_vel );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, YSIZE,  YVEL+4*FSZ, s );

   sprintf( s, "%2.0f", max_vel );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, gright+F4, YVEL+4*FSZ, s );
   
   cx1 =  gscale( angle_vel, min_vel, max_vel ) ;
   cx2 =  gscale(         0, min_vel, max_vel ) ;
   RectE( im, gleft, YVEL+2*FSZ, gright, YVEL+4*FSZ, white );
   RectF( im,   cx1, YVEL+2*FSZ,    cx2, YVEL+4*FSZ, white );

   //  angle gauge -------------------------------

   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0,
                              YSIZE, YANG+2*FSZ-F4, "angle" );
   gdImageStringFT( im, NULL, white, FNT, F2, 0.0,
                              YSIZE, YANG+5*FSZ, "degrees" );

   sprintf( s, "%2.0f", angle_min );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, YSIZE,  YANG+4*FSZ, s );

   sprintf( s, "%2.0f", angle_max );
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, gright+F4, YANG+4*FSZ, s );
   
   cx1 =  gscale(     angle, angle_min, angle_max ) ;
   cx2 =  gscale(         0, angle_min, angle_max ) ;
   RectE( im, gleft, YANG+2*FSZ, gright, YANG+4*FSZ, white );
   RectF( im,   cx1, YANG+2*FSZ,    cx2, YANG+4*FSZ, white );
}



//----------------------------------------
// this draws dots on the upper left edge, useful for scaling sizes in
// openoffice
void lineofdots() {
    int ctr;
    for( ctr=0; ctr < XSIZE; ctr +=4 ) {
       gdImageSetPixel( im, ctr, 0, gray );
       gdImageSetPixel( im, ctr, 1, gray );
       gdImageSetPixel( im, 0, ctr, gray );
       gdImageSetPixel( im, 1, ctr, gray );
    }
}


// SERVERSAT drawing tools
//----------------------------------------
// serversat globals

int       thx, thy ;
int       x,  y                        ;
gdPoint   rbox[4]                      ;
int       light[RIGHT][YSIZE]          ; 
double    xcent, ycent                 ;
int       ly[10]                       ; // levels
int       back[10]                     ; // backside light level

//----------------------------------------
int  y2x( int yy ) {
   return (int)( xcent-angt*(((double)yy)-ycent) ) ;
}

//----------------------------------------
int  f2y( double yy ) {
   return (int)(ycent-angc*yy*scale) ;
}

//---------------------------------------
void chunk( int ya, int yb, int color ) {
   if( ya != yb ) {
      int xa    = y2x(ya)  ;
      int xb    = y2x(yb)  ;
  
      rbox[0].x = xa - thx ;
      rbox[0].y = ya - thy ;
      rbox[1].x = xa + thx ;
      rbox[1].y = ya + thy ;
      rbox[2].x = xb + thx ;
      rbox[2].y = yb + thy ;
      rbox[3].x = xb - thx ;
      rbox[3].y = yb - thy ;
      gdImageFilledPolygon( im, rbox, 4, color ) ;
}  }

//---------------------------------------
//
gdPoint arotate( double rang, int xr, int yr, double rx, double ry ) {
   gdPoint ret;
   ret.x = xr + (int)( rx*cos(rang) - ry*sin(rang) ) ;
   ret.y = yr + (int)( ry*cos(rang) + rx*sin(rang) ) ;

   return ret ;
}

//---------------------------------------
//     x               y          angle        size

void arrow ( int acolor, int xa, int ya, double  aang, double asize ) {
   gdPoint  apoly[7] ;
   if( asize > thick ) {   //  point away
      apoly[0]=arotate( aang, xa, ya, 3.0*thick + asize,  0.0       );
      apoly[1]=arotate( aang, xa, ya, 2.0*thick + asize, -1.0*thick );
      apoly[2]=arotate( aang, xa, ya, 2.0*thick + asize, -0.5*thick );
      apoly[3]=arotate( aang, xa, ya, 2.0*thick        , -0.5*thick );
      apoly[4]=arotate( aang, xa, ya, 2.0*thick        ,  0.5*thick );
      apoly[5]=arotate( aang, xa, ya, 2.0*thick + asize,  0.5*thick );
      apoly[6]=arotate( aang, xa, ya, 2.0*thick + asize,  1.0*thick );
      gdImageFilledPolygon( im, apoly, 7, acolor );
   }
   else if( asize < -thick ) {   //  point towards
      apoly[0]=arotate( aang, xa, ya, 2.0*thick        ,  0.0       );
      apoly[1]=arotate( aang, xa, ya, 3.0*thick        , -1.0*thick );
      apoly[2]=arotate( aang, xa, ya, 3.0*thick        , -0.5*thick );
      apoly[3]=arotate( aang, xa, ya, 3.0*thick - asize, -0.5*thick );
      apoly[4]=arotate( aang, xa, ya, 3.0*thick - asize,  0.5*thick );
      apoly[5]=arotate( aang, xa, ya, 3.0*thick        ,  0.5*thick );
      apoly[6]=arotate( aang, xa, ya, 3.0*thick        ,  1.0*thick );
      gdImageFilledPolygon( im, apoly, 7, acolor );
   }
}

//----------------------------------------
//  scale is Y/2
//  drawn in RIGHT*YSIZE panel

void drawserversat( double t1, double t2 ) {

   int     cline   = black               ; // line color

   scale     = 0.5*SCALE*((double)YSIZE) ; // distance from center to top
   ycent     = 0.55*((double)YSIZE)      ;
   xcent     = 0.5*((double)RIGHT)       ;

   int     xline                         ; // serversat centerline x val
   int     i                             ;
   int     yellow[256]                   ;
   int     xx1, xx2                      ;
   int     yy1, yy2                      ;
   int     yyt, yyb                      ;
   int     xxl, xxr                      ;
   int     xxa                           ;
   int     xy                            ;

   // fill color array  -------

   yellow[0] = gdImageColorAllocate( im,  0,     0,   0 );
   yellow[1] = gdImageColorAllocate( im,  224, 192,   0 );
   yellow[2] = gdImageColorAllocate( im,  255, 224,   0 );
   yellow[3] = gdImageColorAllocate( im,  255, 255, 127 );
   yellow[4] = gdImageColorAllocate( im,  255, 255, 192 );

   // establish l levels ---------------------------
                                           // The below are "pixel" lines
   ly[0]=0                               ; // top of drawing
   back[0]=1 ;
   ly[1]=f2y(  1.0                    )  ; // top of top thruster
   back[1]=1 ;
   ly[2]=f2y(  0.75 +0.25*t2          )  ; // top of top thruster reflective
   back[2]=0 ;
   ly[3]=f2y(  0.75 -0.25*t2          )  ; // bot of top thruster reflective
   back[3]=1 ;
   ly[4]=f2y(  0.5                    )  ; // top of solar cell
   back[4]=0 ;
   ly[5]=f2y( -0.5                    )  ; // bot of solar cell
   back[5]=1 ;
   ly[6]=f2y( -0.75 +0.25*t1          )  ; // top of bot thruster reflective
   back[6]=0 ;
   ly[7]=f2y( -0.75 -0.25*t1          )  ; // bot of bot thruster reflective
   back[7]=1 ;
   ly[8]=f2y( -1.0                    )  ; // bot of bot thruster
   back[8]=1 ;
   ly[9]=YSIZE                           ; // bottom of drawing
   back[9]=1 ;

// printf( "%3.0f%4d%4d%4d%4d%4d%4d%4d%4d%4d%4d\n",
//          angle, thx, thy, ly[1],ly[2],ly[3],ly[4],ly[5],ly[6],ly[7],ly[8] );

   // draw light      ------------------------------

   // work down through the levels
   for( i = 0; i < 9 ; i++ ) { 
      yy1 = ly[ i   ] ;
      yy2 = ly[ i+1 ] ;

      if( yy2 > yy1 ) {
         if(  back[ i ] == 1 ) {
            for( y=yy1 ; y<yy2 ; y++ ) {
               for( x=0 ; x<RIGHT ; x++ ) {
                  light[x][y] = 1 ;
         }  }  }
         else {
            for( y=yy1 ; y<yy2 ; y++ ) {
               xy = y2x( y ) ;
               if( 0 < xy ) { 
                  for( x=0 ; x < xy ; x++ ) {
                     light[x][y] = 1 ;
               }  }
               if( xy < RIGHT ) { 
                  for( x=xy ; x < RIGHT ; x++ ) {
                     light[x][y] = 0    ;
               }  }
   }  }  }  }
      
   // add off axis reflections ---------------------

   for( i = 2 ; i < 9 ; i += 4 ) {   // 2 and 6, cheesy
      if( ly[i] != ly[i+1] ) {       // if beamwidth greater than 0

         yy1 = ly[i]    ;
         yy2 = ly[i+1]  ;
         xx1 = y2x(yy1) ;
         xx2 = y2x(yy2) ;

         //          | bottom |  front |   top  |  back  |
         // an2    -135      -45      +45     +135     -135
         // an2a          -        +        +        -
         // an2b          +        +        -        -

         if( an2a > 0.0 ) {          // front or top
            if( an2b > 0.0 ) {  

            // front -----------------------------------
               if( xx1 > xx2 ) {     // front above
                  xxa = xx1 ;
               } else {              // front below
                  xxa = xx2 ;
               }

               for( x=0 ; x < xxa ; x++ ) {
                  yyt = yy1 - (int)( (an2s/an2c)*(xx1-x) ) ;
                  if( yyt < 0 ) { yyt = 0 ; }
                  yyb = yy2 - (int)( (an2s/an2c)*(xx2-x) ) ;
                  if( yyb > YSIZE-1 ) { yyb = YSIZE-1 ; }
                  if( ( yyb > 0 ) && ( yyt < YSIZE-1) ) { 
                     for( y=yyt ; y <= yyb ; y++ ) {
                        if( x < y2x( y ) ) {
                           light[x][y]++ ;
            }  }  }  }  }  // end of front

            else {                   

            // top   -----------------------------------
               for( y=0 ; y < yy2 ; y++ ) {         // yy2 is always back
                  xxl = xx2 - (int)( (an2c/an2s)*(yy2-y) ) ;
                  if( xxl < 0 ) { xxl = 0 ; }
                  xxr = xx1 - (int)( (an2c/an2s)*(yy1-y) ) ;
                  if( xxr > RIGHT ) { xxr = RIGHT ; }
                  if( (xxr>=0) && (xxl<=RIGHT) )  { 
                     for( x=xxl; x <= xxr ; x++ ) {
                        if( x < y2x( y ) ) {
                           light[x][y]++ ;
         }  }  }  }  }  }  // end of top

         else {                      // back or bottom
            if( an2b > 0.0 ) {       // bottom

            // top   -----------------------------------
               for( y=yy1 ; y < YSIZE ; y++ ) {     // yy1 is always back
                  xxl = xx1 - (int)( (an2c/an2s)*(yy1-y) ) ;
                  if( xxl < 0 ) { xxl = 0 ; }
                  xxr = xx2 - (int)( (an2c/an2s)*(yy2-y) ) ;
                  if( xxr > RIGHT ) { xxr = RIGHT ; }
                  if( (xxr>=0) && (xxl<=RIGHT) )  { 
                     for( x=xxl; x <= xxr ; x++ ) {
                        if( x < y2x( y ) ) {
                           light[x][y]++ ;
            }  }  }  }  }  // end of bottom

            else {                   // back
               // MORE LATER BACK
   }  }  }  }

   // add front light -----------------------------
   for( x=0 ; x<10 ; x++ ) {
      for( y=0 ; y<YSIZE ; y++ ) {
         light[x][y] += 2 ;
   }  }

   // emit light -----------------------------------
   for( x=0 ; x<RIGHT ; x++ ) {
      for( y=0 ; y<YSIZE ; y++ ) {
         gdImageSetPixel( im, x, y, yellow[light[x][y]] ) ;
   }  }

   // draw serversat chunks ------------------------------

   thick = THICK ;
   thx=(int)(0.5*thick*angc)      ;
   thy=(int)(0.5*thick*angs)      ;

   chunk( ly[2], ly[3], white )   ; // draw upper reflector middle
   chunk( ly[4], ly[5], dgray )   ; // draw solar cell
   chunk( ly[6], ly[7], white )   ; // draw lower reflector middle

   // draw serversat lines -------------------------------

   rbox[0].x = y2x(ly[1])-thx ;
   rbox[0].y = ly[1]     -thy ;
   rbox[1].x = y2x(ly[1])+thx ;
   rbox[1].y = ly[1]     +thy ;
   rbox[2].x = y2x(ly[8])+thx ;
   rbox[2].y = ly[8]     +thy ;
   rbox[3].x = y2x(ly[8])-thx ;
   rbox[3].y = ly[8]     -thy ;
   gdImagePolygon( im, rbox, 4, cline );

   // draw light and force arrows
   //     x               y          angle        size
   arrow( white, y2x(f2y( 0.73)), f2y( 0.73),  pi,         -0.3*scale*t2 ) ; 
   arrow( white, y2x(f2y( 0.77)), f2y( 0.77),  pi+an2,      0.3*scale*t2 ) ;
   arrow( white, y2x(f2y( 0.75)), f2y( 0.75),  ang,    0.6*angc*scale*t2 ) ;

   arrow( white, y2x(f2y(-0.77)), f2y(-0.77),  pi,         -0.3*scale*t1 ) ; 
   arrow( white, y2x(f2y(-0.73)), f2y(-0.73),  pi+an2,      0.3*scale*t1 ) ;
   arrow( white, y2x(f2y(-0.75)), f2y(-0.75),  ang,    0.6*angc*scale*t1 ) ;

/*
   // draw fixed pivot arrow
   
   thick = 2*FSZ ;
   yy1 = f2y(0.0) ;
   xx1 = y2x(yy1) - (int)( 2.0*thick ) ;
   arrow( dgray,  xx1, yy1,  0,   -8.0*FSZ ) ;
   yy1 += (int)( 0.5* FSZ );
   xx1  = y2x(yy1) + (int)( 0.5*thick )  ;
   gdImageStringFT( im, NULL, white, FNT, FSZ, 0.0, xx1, yy1, "fixed pivot" );
   thick = THICK ;
*/

}

//=====================================================================

int main () {
   FILE   *pngout;
   char   comment[80]     ;
   char   dirname[80]     ;
   char   framename[80]   ;
   double timem           ;            // time in minutes
   double t1p, t2p        ;
   double t1, t2          ;            // thrust value for each cell
   int    frame           ;            // frame count
   int    npics           ;            // last frame count
   int    i               ; 

   npics  = (10*(TEND-TSTART)+1)/(10*TDELTA) ;  // a bit extra for rounding err
   pi4    = atan(1.0)          ; // pi/4 , 45 degrees
   pi     = 4.0*pi4            ; // pi
   
   // makestars();

   sprintf( rmmkdir, "rm -rf %sdir ; mkdir %sdir", NAME, NAME );

   // make directory
   system( rmmkdir );
 
   for(  frame = 0 ; frame <= npics ; frame++ ) {

      printf( "%04d/%04d\r", frame, npics );
      fflush( stdout );

      im    = gdImageCreateTrueColor(XSIZE, YSIZE);
      white = gdImageColorAllocate (im, 255, 255, 255);
      black = gdImageColorAllocate (im,   0,   0,   0);
      gray  = gdImageColorAllocate (im, 127, 127, 127);
      dgray = gdImageColorAllocate (im,  64,  64,  64);
      lred  = gdImageColorAllocate (im, 255, 208, 208);
      lgreen =gdImageColorAllocate (im, 208, 255, 208);
      trans = gdImageColorAllocate (im,   1,   1,   1);

      thick = THICK ;

      time = TDELTA*( (double) frame ) + TSTART ;
      sprintf( framename, "%sdir/a%04d.png", NAME , frame );
      pngout = fopen ( framename , "wb");

      //  start out stopped, display for a while
      if(   time <=   0.0 ) {
        angle_acc =  0.0 ;
        angle_vel =  0.0 ;
        angle     =  0.0 ;
        t1p       =  0.5 ;
        t2p       =  0.5 ;
        t1        =  0.5 ;
        t2        =  0.5 ;
        spdcolor  =  white ;
        strcpy( comment, "stopped" );
      }

      //  accelerate between   0 seconds and TIME1 seconds 
      else if( time < 1.0*TIME1 ) {
        timem     = time/ MINUTE ;
        angle_acc =  ACC ;
        angle_vel =  ACC * timem ;
        angle     =  ACC * 0.5 * timem * timem ;
        t1p       =  0.05 ;
        t2p       =  0.95 ;
        spdcolor  =  lgreen ;
        strcpy( comment, "speed up" );
      }

      //  decelerate between TIME1 seconds and 2*TIME1 seconds 
      else if( time < 2.0*TIME1 ) {
         timem     = (2.0*TIME1-time)/ MINUTE ;
         angle_acc = -ACC ;
         angle_vel =  ACC * timem ;
         angle     =  angle_max - (0.5*ACC*timem*timem );
         t1p       =  0.95 ;
         t2p       =  0.05 ;
         spdcolor  =  lred ;
         strcpy( comment, "slow down" );
      }
      //  end up stopped, display for a while
      else {
         angle_acc = 0.0 ;
         angle_vel = 0.0 ;
         angle     =  angle_max;
         t1p       =  0.5 ;
         t2p       =  0.5 ;
         spdcolor  =  white ;
         strcpy( comment, "stopped" );
      }

      t1 += FILTER*(t1p-t1 );
      t2 += FILTER*(t2p-t2 );

      ang    = DEG2RAD(  angle )  ; //
      angc   = cos(ang)           ; //
      angs   = sin(ang)           ; //
      angt   = tan(ang)           ; //
      an2    = 2.0*ang            ; // double angle for reflection

      //          | bottom |  front |   top  |  back  |
      // an2    -135      -45      +45     +135     -135
      // an2a          -        +        +        -
      // an2b          +        +        -        -

      an2c   = cos(an2)           ; // 
      an2s   = sin(an2)           ; //
      an2a   = cos(an2-pi4)       ; // +if  -45 < a2 < 135
      an2b   = cos(an2+pi4)       ; // +if -135 < a2 <  45

      // draw serversat -----------------------------------------------
   
      drawserversat( t1, t2 );

      // draw panel and other things ---------------------------------

      lineofdots();                  // calibration strips
      drawgauges();                  // draw gauges

      // status text -------------------------------------------------
      
      /*
      gdImageStringFT( im, NULL, spdcolor,FNT, 2*FSZ, 0.0,
                           FSZ, (int)(2.5*FSZ), comment );
      */
      
      // output frame to directory in PNG format

      gdImagePngEx( im, pngout, 1 );
      gdImageDestroy (im);
      fclose (pngout);
   }

   sprintf( png2swf, "png2swf -o %s.swf -r%3d %sdir/*.png", NAME, rate, NAME );
   system(  png2swf );
}
