/*----------------------------------------*/
/*                                        */
/*     Filename:  hp_math.include         */
/*                                        */
/*     Contains C source code for         */
/*     routines ProcessFloat,             */
/*     ProcessComplex, ProcessInteger,    */
/*     matherr, StoreRegister,            */
/*     RecallRegister, SelectStoreReg     */
/*     and CheckStack.                    */
/*                                        */
/*----------------------------------------*/
ProcessFloat()
{
  double temp;
                                    
  if ( inkey == 52 ) {
    if ( entry_in_progress == FALSE ) {
         t = z;
         z = y;
         y = x;                 }
    return(0);       }


  /* The binary operators that remove y and replace x */
  if ( (inkey ==  4) || (inkey ==  8) || (inkey ==12) ||
       (inkey == 51) || (inkey == 39) )    {
    if ( inkey ==  4 ) x = y * pow( x, -1.0 );
    if ( inkey ==  8 ) x = y*x;
    if ( inkey == 12 ) x = y-x;
    if ( inkey == 51 ) x = y+x;

    if ( inkey == 39 ) x = pow( y, x );
    if ( error_occurred ) return(0);   /* leave the stack alone */
    y = z;
    z = t;
    t = 0.0;
    return(0);                        }

    /*  Polar to rectangular   */
    if ( inkey == 44 )          {
          temp = y * angfac;
          y = x * sin( temp );
          x = x * cos( temp );
          return(0);  }

    /*  Rectangular to polar   */
    if ( inkey == 45 )          {
          temp = y;
          y = atan2( y, x ) / angfac;
          x = sqrt( x*x + temp*temp );
          return(0);                 }


    if ( inkey == 28 ) x =  sin( x   * angfac );
    if ( inkey == 29 ) x =  cos( x   * angfac );
    if ( inkey == 30 ) x =  tan( x   * angfac );
    if ( inkey == 34 ) x = asin( x ) / angfac  ;
    if ( inkey == 35 ) x = acos( x ) / angfac  ;
    if ( inkey == 36 ) x = atan( x ) / angfac  ;


    /* Use math library functions where possible to get friendly errors */
    if ( inkey == 53 ) x = -x;
    if ( inkey == 33 ) x = pow( x, -1.0 ); 
    if ( inkey == 37 ) x = pow( x,  2.0 );
    if ( inkey == 38 ) x = sqrt( x );
    if ( inkey == 40 ) x = log10( x );
    if ( inkey == 41 ) x = log( x ) / log ( 2.0 );
    if ( inkey == 42 ) x = log( x );
    if ( inkey == 46 ) x = pow( 10.0, x );
    if ( inkey == 47 ) x = pow(  2.0, x );
    if ( inkey == 48 ) x = exp( x );

    if ( inkey == 43 ) {
        t = z;
        z = y;
        y = x;
        x = PI;        }


    /* If we've already encountered an error, nevermind. . . */
    if ( error_occurred ) return(0);

    /* Force over(under)flow if number cannot be displayed   */
    if (  abs(x) > 9.99999999999E99 ) x = pow( 10.0,  10000.0 );
    if ( (abs(x) < 1.0000000000E-99) && (abs(x) > 0.0) )
        x = pow( 10.0, -10000.0 );

}
/*--------------------------------------------------------*/
ProcessComplex()
{
  double temp, temp2, r, theta, lnr, coeff, coshx, sinhx, a,b,c,d, xy;
                                    
  if ( inkey == 52 ) {
    if ( entry_in_progress == FALSE ) {
         t = z;
         z = y;
         y = x;                 }
    return(0);       }

  /*   Note that in all the following computations, x represents the     */
  /*   imaginary part of the argument, while y is the real part.  This   */
  /*   notation is contrary to that used in analysis, where the roles    */
  /*   of x and y are precisely reversed.  This arises since it is more  */
  /*   natural to enter the real part of a number first, which causes    */
  /*   it to end up in the y-register.  Just think backwards.  I know    */
  /*   I certainly do.                                                   */


  /* The binary operators that remove t,z and replace y,x */
  if ( (inkey ==  4) || (inkey ==  8) || (inkey ==12) ||
       (inkey == 51) || (inkey == 39) )    {
    if ( inkey ==  4 ) { temp = pow( y*y + x*x , -1.0 );
                         temp2 = y;
                         y = (y*t     + x*z) * temp;
                         x = (temp2*z - x*t) * temp; }
    if ( inkey ==  8 ) { temp = y;
                         y = y*t - x*z;
                         x = x*t + temp*z; }
    if ( inkey == 12 ) { x = z-x;
                         y = t-y; }
    if ( inkey == 51 ) { x = z+x;
                         y = t+y; }

    /*   y to the x     */
    if ( inkey == 39 ) { r = sqrt(t*t + z*z);
                         theta = atan2( z, t );
                         lnr = log(r);
                         coeff = exp( y*lnr - x*theta );
                         temp = y;
                         y = coeff * cos( x*lnr +    y*theta );
                         x = coeff * sin( x*lnr + temp*theta ); }
    t = 0.0;
    z = 0.0;
    return(0);                        }

    /*  Polar to rectangular   */
    if ( inkey == 44 )          {
          temp = y * angfac;
          y = x * sin( temp );
          x = x * cos( temp );
          return(0);  }

    /*  Rectangular to polar   */
    if ( inkey == 45 )          {
          temp = y;
          y = atan2( y, x ) / angfac;
          x = sqrt( x*x + temp*temp );
          return(0);                 }


    if ( ( inkey > 27) && (inkey < 31) ) {
      coshx = cosh( x * angfac );
      sinhx = sinh( x * angfac );
      if (error_occurred) return(0);

      /*   sine   */
      if ( inkey == 28 ) { x =  cos(y*angfac) * sinhx;
                           y =  sin(y*angfac) * coshx; }
      /*   cosine  */
      if ( inkey == 29 ) { x = -sin(y*angfac) * sinhx;
                           y =  cos(y*angfac) * coshx; }
      /*   tangent  */
      if ( inkey == 30 ) { x =  sin(y*angfac);
                           y =  cos(y*angfac);
                           a =  x * coshx;
                           b =  y * sinhx;
                           c =  y * coshx;
                           d = -x * sinhx;
                           temp = pow( c*c + d*d , -1.0 );
                           if ( error_occurred ) return(0);
                           y = (a*c + b*d) * temp;
                           x = (b*c - a*d) * temp; }   }

    /*  arcsine   */
    if ( inkey == 34 ) { xy = x*y;
                         temp = 1.0 - y*y + x*x;
                         r = sqrt( temp*temp + 4.0*xy*xy );
                         theta = atan2( -2.0*xy, temp );
                         if (theta < 0.0) theta = theta + 2.0 * PI;
                         r = sqrt( r );
                         theta = theta / 2.0;
                         b = r * sin( theta ) + y;
                         a = r * cos( theta ) - x;
                         x = -log( sqrt( a*a + b*b ) ) / angfac;
                         y =  atan2( b, a );
                         if ( y < 0.0 ) y = y + 2.0 * PI;
                         y /= angfac; }

    /*  arccosine  */
    if ( inkey == 35 ) { xy = x*y;
                         temp = 1.0 - y*y + x*x;
                         r = sqrt( temp*temp + 4.0*xy*xy );
                         theta = atan2( -2.0*xy, temp );
                         if (theta < 0.0) theta = theta + 2.0 * PI;
                         r = sqrt( r );
                         theta = theta / 2.0;
                         a =  r * cos( theta ) + x;
                         b = -r * sin( theta ) + y;
                         x = -log( sqrt( a*a + b*b ) ) / angfac;
                         y =  atan2( a, b );
                         if ( y < 0.0 ) y = y + 2.0 * PI;
                         y /= angfac; }

    /*  arctangent  */
    if ( inkey == 36 ) { d = (1.0-x)*(1.0-x) + y*y;
                         a = (-x*x - y*y + 1.0) / d;
                         b = (-2.0 * y) / d;
                         r = sqrt( a*a + b*b );
                         theta = atan2( b, a );
                         y = -theta / 2.0 / angfac;
                         x = log(r) / 2.0 / angfac; }


    if ( inkey == 53 ) { x = -x;
                         y = -y; }
    /*    1/x      */
    if ( inkey == 33 ) { r = pow( x*x + y*y, -1.0 );
                         y =  r*y;
                         x = -r*x; }
    /*    x squared   */
    if ( inkey == 37 ) { r = y;
                         y = y*y - x*x;
                         x = 2.0*x*r; }
    /*    sqrt(x)     */
    if ( inkey == 38 ) { theta = atan2( x, y ) / 2.0;
                         r = sqrt( sqrt( x*x + y*y ) );
                         x = r * sin( theta );
                         y = r * cos( theta );  }

    /*   logarithms   */
    if ((inkey > 39) && (inkey < 43)) {
      theta = atan2( x, y );
      r = sqrt( x*x + y*y );
      if ( inkey == 40 ) a = log( 10.0 );
      if ( inkey == 41 ) a = log(  2.0 );
      if ( inkey == 42 ) a = 1.0;
      x = theta  / a;
      y = log(r) / a;
      }

    /*  exponentials  */
    if ((inkey > 45) && (inkey < 49)) {
      if ( inkey == 46 ) a = log( 10.0 );
      if ( inkey == 47 ) a = log(  2.0 );
      if ( inkey == 48 ) a = 1.0;
      r = exp( y * a );
      y = r * cos( x * a );
      x = r * sin( x * a );
      }

    if ( inkey == 43 ) {
        t = z;
        z = y;
        y = x;
        x = PI;        }


    /* If we've already encountered an error, nevermind. . . */
    if ( error_occurred ) return(0);

    /* Force over(under)flow if number cannot be displayed   */
    if (  abs(x) > 9.99999999999E99 ) x = pow( 10.0,  10000.0 );
    if ( (abs(x) < 1.0000000000E-99) && (abs(x) > 0.0) )
        x = pow( 10.0, -10000.0 );

    /* Force over(under)flow if number cannot be displayed   */
    if (  abs(y) > 9.99999999999E99 ) y = pow( 10.0,  10000.0 );
    if ( (abs(y) < 1.0000000000E-99) && (abs(y) > 0.0) )
        y = pow( 10.0, -10000.0 );
}
/*--------------------------------------------------------*/
ProcessInteger()
{
  x = ix;
  y = iy;
  z = iz;
  t = it;

  ProcessFloat();

  /*  Truncate toward zero after slight rounding:        */
  ix = ( x > 0.0 )  ?  x + TRUNC_MARGIN : x - TRUNC_MARGIN;
  iy = ( y > 0.0 )  ?  y + TRUNC_MARGIN : y - TRUNC_MARGIN;
  iz = ( z > 0.0 )  ?  z + TRUNC_MARGIN : z - TRUNC_MARGIN;
  it = ( t > 0.0 )  ?  t + TRUNC_MARGIN : t - TRUNC_MARGIN;


}
/*-------------------------------------------------------*/
matherr( exc )
struct exception *exc;
{
  if (error_occurred) return(0);
  error_occurred = TRUE;

  DrawDisplay();

  Move( rp, 250, 55 );
  Text( rp, errstring[ exc->type ], 14 );

  Acknowledge:
  WaitPort( hp_window -> UserPort );

  if ( GadgetPoked() )        {
    CloseWindow( hp_window );
    CloseScreen( hp_screen );
    if ( print_on )         {
      fputs("\33#1", printer );    /* restore default setting */
      fclose( printer );    }
    exit();                   }

  if ( message -> Code != SELECTDOWN ) goto Acknowledge;

  DrawDisplay();

  /*  Zero out appropriate stack elements  */
  if ( base < BINARY )   {
    x = 0.0;
    if ( base == COMPLEX ) y = 0.0;
    DisplayFloatXY();    }
  else                   {
    ix = 0;
    DisplayIntXY();      }

}
/*------------------------------------------------------*/
StoreRegister()
{
  if ( entry_in_progress ) {
       entry_in_progress = FALSE;
       PushX();            }

  Move( rp, start_display[ base ], 59 );
  Text( rp, "                                ", display_length[ base ] );
  Move( rp, start_display[ base ], 59 );
  Text( rp, "  STO  _ ", 9 );

  SelectStoreReg();

  if ( base < BINARY )  {
    registers[reg] = x;
    DisplayFloatXY();   }
  else                  {
    iregisters[reg] = ix;
    DisplayIntXY();     }

}
/*-----------------------------------------------------------*/
RecallRegister()
{
  if ( entry_in_progress ) {
       entry_in_progress = FALSE;
       PushX();            }

  Move( rp, start_display[ base ], 59 );
  Text( rp, "                                ", display_length[ base ] );
  Move( rp, start_display[ base ], 59 );
  Text( rp, "  RCL  _ ", 9 );

  SelectStoreReg();

  if ( base < BINARY )  {
    t = z;
    z = y;
    y = x;
    x = registers[reg];
    DisplayFloatXY();   }
  else                  {
    it = iz;
    iz = iy;
    iy = ix;
    ix = iregisters[reg];
    DisplayIntXY();     }

}
/*-----------------------------------------------------------*/
SelectStoreReg()
{
  int k;

  Sleep:
  WaitPort( hp_window -> UserPort );

  if ( GadgetPoked() )        {
    CloseWindow( hp_window );
    CloseScreen( hp_screen );
    if ( print_on )       {
      fputs("\33#1", printer );    /* restore default setting */
      fclose( printer );  }
    exit();                   }

  if ( message -> Code != SELECTDOWN )  goto Sleep;

  regcode = KeyCode();

  if ( (regcode > 18) && (regcode < 49)  )  goto Sleep;
  if (  regcode > 49)                       goto Sleep;
  if ( (regcode == 0) || (regcode ==  4) ||
       (regcode == 8) || (regcode == 12  )) goto Sleep;

  
  Move( rp, start_display[ base ] + 39, 59 );
  Text( rp, opcode[ regcode ], 5 );

  for ( k = 1; k < 7000; k++) xmin = 5;

  if ( regcode ==  1 ) reg = 7;
  if ( regcode ==  2 ) reg = 8;
  if ( regcode ==  3 ) reg = 9;
  if ( regcode ==  5 ) reg = 4;
  if ( regcode ==  6 ) reg = 5;
  if ( regcode ==  7 ) reg = 6;
  if ( regcode ==  9 ) reg = 1;
  if ( regcode == 10 ) reg = 2;
  if ( regcode == 11 ) reg = 3;
  if ( (regcode > 12) && (regcode <19) ) reg = regcode - 3;
  if ( regcode == 49 ) reg = 0;

}
/*-----------------------------------------------------*/
CheckStack()
{
  double t_float;
  int    t_int;

  if ( inkey < 56 ) return(0);

  /*     RCLz     */
  if ( inkey == 56 ) {  it = iz;
                        iz = iy;
                        iy = ix;
                        ix = it;   /* was iz a moment ago */
                        t  =  z;
                        z  =  y;
                        y  =  x;
                        x  =  t;   }

  /*     RCLt     */
  if ( inkey == 57 ) {  t_int =  it;
                        it    =  iz;
                        iz    =  iy;
                        iy    =  ix;
                        ix    =  t_int;
                        t_float = t;
                        t       = z;
                        z       = y;
                        y       = x;
                        x       = t_float; }

  /*      LastX            */
  if ( inkey == 58 )     {
      tempbase = base;
      base = lastbase;
      PushX();
      base = tempbase;   }

  /*      x <> y           */
  if ( inkey == 59 ) {  t_int   =      ix;
                        ix      =      iy;
                        iy      =   t_int;
                        t_float =       x;
                        x       =       y;
                        y       = t_float;   }
  

  if ( (inkey == 60) && ( entry_in_progress == FALSE ) ) {
     x  =   y;
     y  =   z;
     z  =   t;
     t  = 0.0;
     ix =  iy;
     iy =  iz;
     iz =  it;
     it =   0;  }

}