tomwalters@0: /*
tomwalters@0:     Copyright (c) Applied Psychology Unit, Medical Research Council. 1988, 1989
tomwalters@0:     ===========================================================================
tomwalters@0: 
tomwalters@0:     Permission to use, copy, modify, and distribute this software without fee 
tomwalters@0:     is hereby granted for research purposes, provided that this copyright 
tomwalters@0:     notice appears in all copies and in all supporting documentation, and that 
tomwalters@0:     the software is not redistributed for any fee (except for a nominal shipping 
tomwalters@0:     charge). Anyone wanting to incorporate all or part of this software in a
tomwalters@0:     commercial product must obtain a license from the Medical Research Council.
tomwalters@0: 
tomwalters@0:     The MRC makes no representations about the suitability of this 
tomwalters@0:     software for any purpose.  It is provided "as is" without express or implied 
tomwalters@0:     warranty.
tomwalters@0:  
tomwalters@0:     THE MRC DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING 
tomwalters@0:     ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS, IN NO EVENT SHALL THE
tomwalters@0:     A.P.U. BE LIABLE FOR ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY 
tomwalters@0:     DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN 
tomwalters@0:     AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF 
tomwalters@0:     OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
tomwalters@0: */
tomwalters@0: 
tomwalters@0: #include <stdio.h>
tomwalters@0: #include <math.h>
tomwalters@0: 
tomwalters@0: #include "stitch.h"
tomwalters@0: 
tomwalters@0: #define OUTPUT_SCALE 360
tomwalters@0: 
tomwalters@0: #ifdef TEST
tomwalters@0: 
tomwalters@0: #define TwoPi (3.1415926535*2)
tomwalters@0: 
tomwalters@0: 
tomwalters@0: main( argc, argv )
tomwalters@0: int   argc ;
tomwalters@0: char *argv[] ;
tomwalters@0: {
tomwalters@0:     double deg ;
tomwalters@0: 
tomwalters@0:     printf( "hello\n" );
tomwalters@0: 
tomwalters@0:     for( deg = 0 ; deg < TwoPi/2 ; deg += TwoPi/90001. )
tomwalters@0: 	if( (int) ( atan2( sin( deg ), cos( deg ) )/TwoPi*OUTPUT_SCALE + 0.5 ) != iatan2( (int) ( sin(deg)*1000 ), (int) ( cos(deg)*1000 ) ) )
tomwalters@0: 	    printf( " %f %d\n", atan2( sin( deg ), cos( deg ) )/TwoPi*OUTPUT_SCALE, iatan2( (int) ( sin(deg)*1000 ), (int) ( cos(deg)*1000 ) ) ) ;
tomwalters@0: 
tomwalters@0:     for( ; deg < TwoPi ; deg += TwoPi/90001. )
tomwalters@0: 	if( (int) ( atan2( sin( deg ), cos( deg ) )/TwoPi*OUTPUT_SCALE - 0.5 ) != iatan2( (int) ( sin(deg)*1000 ), (int) ( cos(deg)*1000 ) ) )
tomwalters@0: 	    printf( " %f %d\n", atan2( sin( deg ), cos( deg ) )/TwoPi*OUTPUT_SCALE, iatan2( (int) ( sin(deg)*1000 ), (int) ( cos(deg)*1000 ) ) ) ;
tomwalters@0: 
tomwalters@0: 
tomwalters@0: }
tomwalters@0: 
tomwalters@0: #endif
tomwalters@0: 
tomwalters@0: iatan2( y, x )
tomwalters@0: int     y, x ;
tomwalters@0: {
tomwalters@0:     static int table_bits = 14 ;
tomwalters@0:     static unsigned table_size ;
tomwalters@0:     static short *atan_table = ( short * ) 0 ;
tomwalters@0:     double two_phi ;
tomwalters@0:     int i ;
tomwalters@0: 
tomwalters@0:     if( atan_table == ( short * ) 0 ) {
tomwalters@0: 
tomwalters@0: 	table_size = 1 << table_bits ;
tomwalters@0: 
tomwalters@0: 	atan_table = ( short * ) stitch_malloc( sizeof ( *atan_table ) * ( table_size + 1 ), "Making atan lookup table" ) ;
tomwalters@0: 
tomwalters@0: 	two_phi = atan( 1. ) * 8. ;
tomwalters@0: 
tomwalters@0: 	for( i=0 ; i <= table_size ; i++ )
tomwalters@0: 	    atan_table[ i ] = atan( ( double ) i / ( table_size ) ) / two_phi * OUTPUT_SCALE + 0.5 ;
tomwalters@0:     }
tomwalters@0: 
tomwalters@0:     if( x != 0 || y != 0 )
tomwalters@0: 	if( y > 0 )
tomwalters@0: 	    if( x > 0 )
tomwalters@0: 		if(  x >  y )
tomwalters@0: 		    return(                          atan_table[ ( (  y << table_bits + ( table_size >> 1 ) ) ) /  x ] ) ;
tomwalters@0: 		else
tomwalters@0: 		    return(  ( OUTPUT_SCALE >> 2 ) - atan_table[ ( (  x << table_bits + ( table_size >> 1 ) ) ) /  y ] ) ;
tomwalters@0: 	    else
tomwalters@0: 		if( -x <  y )
tomwalters@0: 		    return(  ( OUTPUT_SCALE >> 2 ) + atan_table[ ( ( -x << table_bits + ( table_size >> 1 ) ) ) /  y ] ) ;
tomwalters@0: 		else
tomwalters@0: 		    return(  ( OUTPUT_SCALE >> 1 ) - atan_table[ ( (  y << table_bits + ( table_size >> 1 ) ) ) / -x ] ) ;
tomwalters@0: 	else
tomwalters@0: 	    if( x < 0 )
tomwalters@0: 		if( -x > -y )
tomwalters@0: 		    return( -( OUTPUT_SCALE >> 1 ) + atan_table[ ( ( -y << table_bits + ( table_size >> 1 ) ) ) / -x ] ) ;
tomwalters@0: 		else
tomwalters@0: 		    return( -( OUTPUT_SCALE >> 2 ) - atan_table[ ( ( -x << table_bits + ( table_size >> 1 ) ) ) / -y ] ) ;
tomwalters@0: 	    else
tomwalters@0: 		if(  x < -y )
tomwalters@0: 		    return( -( OUTPUT_SCALE >> 2 ) + atan_table[ ( (  x << table_bits + ( table_size >> 1 ) ) ) / -y ] ) ;
tomwalters@0: 		else
tomwalters@0: 		    return(                        - atan_table[ ( ( -y << table_bits + ( table_size >> 1 ) ) ) /  x ] ) ;
tomwalters@0:     else
tomwalters@0: 	return( 99 ) ;
tomwalters@0: }