/* 
 * PlotView.m -- Implementation file for the PlotView class 
 *
 * You may freely copy, distribute, and reuse the code in this example.
 * NeXT disclaims any warranty of any kind, expressed or implied, as to its
 * fitness for any particular use.
 *
 * This is an implementation of the PlotView class to demonstrate the effects 
 * and use of simulated annealing, in particular the Traveling Salesman Problem
 * (TSP). Implemented by Nolan R. Davis and Joseph Jeffery in the Acoustics 
 * Division of the Naval Research Laboratory.
 */

#import "PlotView.h"
#import "drawingFuncs.h"
#import <objc/Storage.h>
#import <appkit/NXCursor.h>
#import <appkit/Application.h>
#import <appkit/Window.h>
#import <appkit/FormCell.h>
#import <appkit/Control.h>
#import <appkit/Pasteboard.h>
#import <NXCType.h> 
#import <math.h> 
#import <dpsclient/psops.h>
#import <dpsclient/wraps.h>
#import <stdlib.h> 
#import <string.h> 
#define MAXNUMLENGTH 50
#define DEFAULTRADIUS 4.0
int route = 0;
float mu = 1.0;
  long		ialpha;
  long		ibeta;
  long		ir;
int ttt=0;
/* Subroutines used in the annealing, found at the end. */
extern double get_distance( int k , int l );
extern double min_value ( double a, double b );
extern double get_uniform_random( );
extern float getNumber(NXStream *stream);
extern int getSeparator(NXStream *stream);

/* Structures used to store the coordinates of the cities and the paths between them. */
  struct	 {
	double	x;
	double	y;
}city[30];

 struct	{
	int 	city;
	double	length;
}path[30],bestpath[30];


@implementation PlotView

- initFrame:(const NXRect *)frameRect

/* Initializes the new PlotView object. First, an initFrame: message is sent to  super to initialize PlotView as a View.  Next, the PlotView sets its own state -- that it is opaque and that the origin of its coordinate system lies in  the center of its area. It then creates and initializes its associated objects, a Storage object, an NXCursor, and a Cell.  Finally, it loads into the Window Server some PostScript procedures that it will use in drawing itself. 
 */

{
	NXPoint	spot;
	
	[super initFrame:frameRect];
	[self setOpaque:YES];
	[self setRadius:DEFAULTRADIUS];
	[self translate:floor(frame.size.width/2) :floor(frame.size.height/2)];
	[self scale:1 :1];



	points = [[Storage alloc] initCount:0 elementSize:sizeof(NXPoint) description:"{ff}"];
	crossCursor = [NXCursor newFromImage:[NXImage newFromSection:"cross.tiff"]];
	spot.x = spot.y = 7.0;
	[crossCursor setHotSpot:&spot];
	/* 
	 * Normally, the next two message expressions would be combined, as:
	 * 		readOut = [[Cell alloc] initTextCell:""];
	 * but are here separated for clarity.
	 */
	/*readOut = [[FormCell alloc] initTextCell:""];*/
	loadPSProcedures();
	return self;
}

- setTextCell:aTextField
{
readOut=aTextField;	

return self;
}
- setDelegate:anObject
/*
 * Sets the PlotView's delegate instance variable to the supplied object.
 */
{
	delegate = anObject;
	return self;
}


- drawSelf:(const NXRect *)rects :(int)rectCount
/*
 * Draws the PlotView's background and axes.  If there are any points,
 * these are drawn too.
 *
 * (Note:  For simplicity, although PlotView only repaints the background 
 * of the update region, it redraws the entire axes.  For better performance, 
 * only those parts of the axes that fall within the update region should 
 * be redrawn.)   
 */
{
	unsigned int	 i,j;
	NXPoint		*aPoint, *bPoint;
	int c;
	if (rects == NULL) return self;

	/* If we're printing, we need to load drawing procedures to the printing context.  This simple test is OK for a View that prints a single page.  If there were multiple pages, the procedures would be loaded multiple times.  Although this is inefficient, it causes no error.  Views that print multiple pages should load procedures such as these in an endPrologue method (see the View spec sheet). */

	if (NXDrawingStatus != NX_DRAWING) 
		loadPSProcedures();

	/* paint visible area white then draw axes */
	PSsetgray(NX_WHITE);
	NXRectFill(&rects[0]);
	PSsetgray(NX_DKGRAY);
	drawAxes(bounds.origin.x, bounds.origin.y, bounds.size.width, bounds.size.height);
        c=0;
	
	/* now take each point and draw it */
	PSsetgray(NX_BLACK);
	i = [points count];
       j=i;
	while (i--) {
		aPoint = (NXPoint *)[points elementAt:i];
		drawCircle(aPoint->x, aPoint->y, radius);
		
	/* Now draw the connecting lines */	
		if(route==1){
		  if(c==1){
		    bPoint = (NXPoint *)[points elementAt:i+1];
		    join(aPoint->x, aPoint->y,bPoint->x, bPoint->y);}
		    else{bPoint = (NXPoint *)[points elementAt:0];join(aPoint->x, aPoint->y,bPoint->x, bPoint->y);}c=1;}
        } 
	 
	return self;
}

- DoRoute:sender
{
/* Called when the Find & Show Route button is pressed. It does the main part of the annealing  and then calls plot self to display it*/
  unsigned      int	 i;
  NXPoint		*aPoint, bPoint;
  int           NCITY,itmp,n_iterations,n_no_change,change,z,k3,trials;
  static	int	j;
  static	int	k1;
  static	int	k2;
  static	double	e;
  static	double	delta_e;
  static	double	best_e;
  static	double	temperature;
  static	double	e_baseline;
  static	double	r;
  static	double 	trial_length_1_in;
  static	double 	trial_length_1_out;
  static	double 	trial_length_2_in;
  static	double 	trial_length_2_out;
  static	double 	length_1_in;
  static	double 	length_2_out;
  static	double 	length_2_in;
  static	double 	length_1_out;
  ialpha = pow(7,5);
  ibeta = pow(2,13) - 1;
  ir = 111;
route=1;

/* Read in the values from the window */
[initTemp takeIntValueFrom:initTemp];
temperature =[initTemp intValue];
[Iterations takeIntValueFrom:Iterations];
trials =[Iterations intValue]+1;
n_iterations = 2;
[freeze takeIntValueFrom:freeze];
change = [freeze intValue];

n_no_change =0;
i = [points count];
NCITY = i;
e=0;
z=0;
ttt=1;
while (i--) {aPoint = (NXPoint *)[points elementAt:i];
     city[i].x=aPoint->x;city[i].y=aPoint->y;
}

for ( j=0; j < NCITY; j++) {
     path[j].city = j ;
     k1 = path[j].city;
     path[j].length = get_distance( j, (j+1)%NCITY );
     e += path[j].length;
     bestpath[j].city= path[j].city;
  } 
best_e=e;e_baseline=e;
while (n_no_change < change) {

  n_iterations++;
	
  for ( k1 = 0; k1 < NCITY; k1++ ) {
	
      do {
         r = get_uniform_random();
         k2 = floor( NCITY * r );
      }
      while (k1==k2);	z++;           
      
      /* Plots out intermediate steps of the annealing if requested. */
      
      if (trials>0){
         if(z%trials==0){ for ( j=0; j <= NCITY-1; j++) {
            k3 = bestpath[j].city;
            bPoint.x =city[k3].x;
            bPoint.y =city[k3].y;
            [points replace:&bPoint at:j];
         }	
        [delegate NukeCities:self];
        for(j=0;j<=NCITY-1;j++) {
            if (delegate && [delegate respondsTo:@selector(plotView:pointDidChange:)])
            [delegate plotView:self pointDidChange:[points elementAt:j]];
            }
        [self plot:self];}
     }

 /* The possible new path lengths */
 
      trial_length_2_in = get_distance ( path[k2].city, path[ (k1-1+NCITY)%NCITY ].city );
      trial_length_1_out = get_distance ( path[k1].city, path[ (k2+1)%NCITY ].city );
      trial_length_2_out = get_distance ( path[k2].city, path[ (k1+1+NCITY)%NCITY ].city );
      trial_length_1_in = get_distance ( path[k1].city, path[ (k2-1+NCITY)%NCITY ].city );
   
   /* The old path lengths*/
   
      length_1_out = get_distance ( path[k1].city, path[ (k1-1+NCITY)%NCITY ].city );
      length_1_in = get_distance ( path[k1].city, path[ (k1+1)%NCITY ].city );
      length_2_out = get_distance ( path[k2].city, path[ (k2-1+NCITY)%NCITY ].city );
      length_2_in = get_distance ( path[k2].city, path[ (k2+1)%NCITY ].city );

       /* Compute energy change */
/* The case where k2 differs from k1 by 1 requires special handling.  */
/* The case where k2 differs from k1 by 2 is more efficient with special handling. */
/* Also, this particular representation depends on which of k1 and k2 is greater.   */
 
      if ( k2 == (k1+1)%NCITY ) {
              delta_e = trial_length_1_out + trial_length_2_in
                                    - length_1_out - length_2_in;
	         if ( get_uniform_random() < exp(-delta_e/temperature) ) {

         /* The state is accepted. Interchange intermediate cities */

        e += delta_e;
        itmp = path[k1].city;
	path[k1].city = path[k2].city;
	path[k2].city = itmp;
    
	if(best_e>e){for ( j=0; j <= NCITY-1; j++) {
        bestpath[j].city= path[j].city;}best_e=e;}
}
}
      else if (k2 == (k1-1+NCITY)%NCITY  )  {
               delta_e = trial_length_1_in + trial_length_2_out
                                    - length_1_in - length_2_out;
	            if ( get_uniform_random() < exp(-delta_e/temperature) ) {

         /* The state is accepted. Interchange intermediate cities */

        e += delta_e;
        itmp = path[k1].city;
	path[k1].city = path[k2].city;
	path[k2].city = itmp;

	if(best_e>e){for ( j=0; j <= NCITY-1; j++) {
        bestpath[j].city= path[j].city;}best_e=e;}
}
}
           else {
              delta_e = trial_length_1_in + trial_length_1_out + trial_length_2_in + trial_length_2_out
                    - length_1_in - length_1_out- length_2_in - length_2_out;
	            if ( get_uniform_random() < exp(-delta_e/temperature) ) {

         /* The state is accepted. Interchange intermediate cities */

        e += delta_e;
	
        itmp = path[k1].city;
	path[k1].city = path[k2].city;
	path[k2].city = itmp;
        
if(best_e>e){for ( j=0; j <= NCITY-1; j++) {
        bestpath[j].city= path[j].city;}best_e=e;}
	

}
}


  }	
    /* Have finished one loop through the cities.	*/
    /* Now reduce the temperature.			*/


    temperature *= min_value( ( (float) n_iterations ) /  ( (float) (n_iterations+1) ) , 
                     2./(1.+exp(mu*(e_baseline/e - 1.) ) ) );


[temp setIntValue:temperature];

    if (e < e_baseline) {
      e_baseline = e;
      n_no_change = 0;
[pathL setIntValue:e_baseline];


    }
    else n_no_change++;

  }
  /* End of while loop.	*/



  /* At this point the system is frozen.  Plot results (best path not final path).*/
printf("Z=%d",z);
  printf( "   Final energy:  %8f\n",best_e);
  for ( j=0; j <= NCITY-1; j++) {
    k3 = bestpath[j].city;
    printf  ( "Final:    %8f %8f\n", city[k3].x, city[k3].y);
    bPoint.x =city[k3].x;
 bPoint.y =city[k3].y;
 [points replace:&bPoint at:j];

  }	
[delegate NukeCities:self];
 for(j=0;j<=NCITY-1;j++) {
        if (delegate && [delegate respondsTo:@selector(plotView:pointDidChange:)])
                    [delegate plotView:self pointDidChange:[points elementAt:j]];
        }
[self plot:self];
route=0;
ttt=0;
return self;
}


- clear:sender
/*
 * Clears the PlotView by emptying its Storage object of all points
 * and then redisplaying the PlotView.  This action is taken only if
 * the PlotView's state requires it.
 */
{
	if (needsClearing) {
		[points empty];
		[self display];
		needsClearing = NO;
	}
	return self;
}

- sizeTo:(NXCoord)width :(NXCoord)height
/*
 * Ensures that whenever the PlotView is resized, the origin of its 
 * coordinate system is repositioned to the center of its area.
 */
{
	[super sizeTo:width :height];
	[self setDrawOrigin:-floor(width/2) : -floor(height/2)]; 
	return self;
}

- registerPoint:(NXPoint *)aPoint
/*
 * Adds a point to the list the PlotView keeps in its Storage object. 
 */
{
	[points addElement:aPoint];
	return self;
}

- mouseDown:(NXEvent *) theEvent
/*
 * Responds to a message the system sends whenever the user presses the mouse
 * button when the cursor is over the PlotView.  The PlotView changes the 
 * cursor to a cross-hairs image and then starts asking for mouse-dragged or mouse-
 * up events.  As it receives mouse-dragged events, the PlotView updates the readOut
 * text Cell with the cursor's coordinates.  If the user releases the mouse
 * button while the cursor is over the PlotView, the PlotView registers the 
 * point and then sends a message to its delegate notifying it of the new 
 * point. If there is a city at the location of the mouseDown its value is 
 * eliminated from the list while the new mouse dragged location is put in its
 * place. Effectually grabbing and moving the point to a new location.
 * 
 */
{
	int	looping = YES, oldMask,distsq,i,X,Y,X1,Y1,z,j;
	NXPoint	aPoint,*cPoint,*bPoint;
	NXRect	plotRect, cellRect;
   	char	buffer[100];


 [crossCursor set];
 [self getBounds:&plotRect];
 NXSetRect(&cellRect, plotRect.origin.x, plotRect.origin.y, 100.0, 20.0);
 oldMask = [window addToEventMask:NX_MOUSEDRAGGEDMASK];
 [self lockFocus];
 aPoint = theEvent->location;
 [self convertPoint:&aPoint fromView:nil];
 distsq=0;
 z=1;
 i = [points count];

 /* See if the new mouse down is at the same location as an old city. */
 
 while (i--) {
    if(z==1){bPoint = (NXPoint *)[points elementAt:i];
         X=bPoint->x;Y=bPoint->y;
         cPoint = &aPoint;
         X1=cPoint->x;Y1=cPoint->y;
         distsq = (X-X1)*(X-X1)+(Y-Y1)*(Y-Y1);
         if(distsq<11){z=0;j=i;}}
 }
 if(z==0){
 
 /* The mouse was over another city, eliminate it and show new city. */
 
	do {
	   aPoint = theEvent->location;
	   [self convertPoint:&aPoint fromView:nil];
	   sprintf(buffer, "(%d,%d)->(%d, %d)", X,Y,(int)aPoint.x, (int)aPoint.y); 
  	   [readOut setStringValue:buffer];
 	   [window flushWindow];
	   if (theEvent->type == NX_MOUSEUP) {
              /*  on mouse-up, register point, inform delegate, and clean up state */
		[readOut setStringValue:""];
  	        /*[readOut drawInside:&cellRect inView:self];*/
		[window flushWindow];
		if (NXPointInRect(&aPoint, &plotRect)) {
			PSsetgray(NX_BLACK);
			drawCircle(aPoint.x, aPoint.y, radius);
			[window flushWindow];
			bPoint = &aPoint;
			[points replace:&aPoint at:j];
			[self display];
                        [delegate NukeCities:self];
	                i = [points count];
	                for(j=0;j<=i-1;j++) {

                             if (delegate && [delegate respondsTo:@selector(plotView:pointDidChange:)])
                            [delegate plotView:self pointDidChange:[points elementAt:j]];
                        }
			needsClearing = YES;
		}
		looping = NO;
  	    }	
   	}while(looping &&(theEvent=[NXApp getNextEvent:NX_MOUSEUPMASK|NX_MOUSEDRAGGEDMASK]));
	[self unlockFocus];
	[window setEventMask:oldMask];
	[NXArrow set];
	return self;
}
else{

/* The mouse was not over another city: treat normally and add the new city to the map and memory. */

	do {
		aPoint = theEvent->location;
		[self convertPoint:&aPoint fromView:nil];
		sprintf(buffer, "(%d, %d)", (int)aPoint.x, (int)aPoint.y); 
  		[readOut setStringValue:buffer];
  		/*[readOut drawInside:&cellRect inView:self];*/
		[window flushWindow];
		if (theEvent->type == NX_MOUSEUP) {
	    /* on mouse-up, register point, inform delegate, and clean up state */
			[readOut setStringValue:""];
  			/*[readOut drawInside:&cellRect inView:self];*/
			[window flushWindow];
			if (NXPointInRect(&aPoint, &plotRect)) {
				PSsetgray(NX_BLACK);
				drawCircle(aPoint.x, aPoint.y, radius);
				[window flushWindow];
				[self registerPoint:&aPoint];
				needsClearing = YES;
				if (delegate && [delegate respondsTo:@selector(plotView:pointDidChange:)])
					[delegate plotView:self pointDidChange:&aPoint];
			}
			looping = NO;
		}	
   	} while(looping &&(theEvent=[NXApp getNextEvent:NX_MOUSEUPMASK|NX_MOUSEDRAGGEDMASK]));
	[self unlockFocus];
	[window setEventMask:oldMask];
	[NXArrow set];
	return self;
}}

- Nuke:sender
{

/* Clears all memory of cities, clears the map and the list.*/

  [delegate NukeCities:sender];
  [self clear:self];
  return self;
}

- setRadius:(float)aFloat
/*
 * Sets the value for the radius of the PlotView's points.
 */
{
	radius = aFloat;
	return self;
}

- setXmax:anObject
{
    xmax = anObject;
    return self;
}
- setXmin:anObject
{
    xmin = anObject;
    return self;
}
- setYmax:anObject
{
    ymax = anObject;
    return self;
}
- setYmin:anObject
{
    ymin = anObject;
    return self;
}
- setFreeze:anObject
{
    freeze = anObject;
    return self;
}
- setInitTemp:anObject
{
    initTemp = anObject;
    return self;
}
- setPathL:anObject
{
    pathL = anObject;
    return self;
}
- setTemp:anObject
{
    temp = anObject;
    return self;
}


- (float)radius
/*
 * Returns the value for the radius of the PlotView's points.
 */
{
	return radius;
}


- read:(NXTypedStream *)stream
/*
 * Unarchives the PlotView.  Initializes four of the PlotView's instance variables
 * to the values stored in the stream. 
 */
{
	[super read:stream];
	delegate = NXReadObject(stream);
	NXReadTypes(stream, "@@@f", &points, &crossCursor, &readOut, &radius);
	return self;
}

- write:(NXTypedStream *)stream
/*
 * Archives the PlotView by writing its important instance variables to the stream. 
 */
{
	[super write:stream];
	NXWriteObjectReference(stream, delegate);
	NXWriteTypes(stream, "@@@f", &points, &crossCursor, &readOut, &radius);
	return self;
}

- awake
/*
 * Finishes initializing a PlotView that has just been unarchived (see the read: method).
 */
{
	loadPSProcedures();
	return [super awake];
}


- (const char *)inspectorName
{
    return "PlotViewInspector";
}

- CoolingOff:sender
{
mu =1.0;
return self;
}
- CoolingOn:sender
{
mu =0.0;
return self;
}


- plot:sender
/*
 * Responds to a plot: message by asking the PlotView's delegate for a stream
 * containing the points to plot.  It then extracts number pairs from the stream
 * and creates NXPoint structures with them.  For each structure it creates, the 
 * PlotView sends itself a registerPoint: message to add the point to its list of 
 * points.  This simple parser ignores input it doesn't understand.
 */
{
	int	 	 i,NCITY,c, sign, retval ,pathlength;
	float	 	 aNum;
	NXPoint	 	 aPoint;
  	NXPoint		*bPoint;
	NXStream	*stream;
	BOOL		 processedLine = NO, gotFirst = NO, gotSecond = NO;
	route=1;
		[self clear:self];
	if (delegate && [delegate respondsTo:@selector(plotView:providePoints:)])
		[delegate plotView:self providePoints:&stream];
	NXSeek(stream, 0, NX_FROMSTART);
	while ((c = NXGetc(stream)) != EOF){
		sign = 1;
		retval = getSeparator(stream);
		if (retval == 1) {
			c = NXGetc(stream);
		} else if (retval == -1) {
			break;
		}
		if (!NXIsDigit(c)) {
			if ((c == '-') && NXIsDigit(c = NXGetc(stream))){
					sign = -1;
			} else {
				while( (c != '\n') && (c != EOF)) {
					c = NXGetc(stream);
					processedLine = YES;
				}
			}
		}
		if (c == EOF)
			break;
		else if (processedLine) {
			processedLine = NO;
			continue;
		}
		aNum = getNumber(stream);
		if (!gotFirst) {
			aPoint.x = sign * aNum;
			gotFirst = YES;
		} else if (!gotSecond) {
			aPoint.y = sign * aNum;
			[self registerPoint:&aPoint];
			gotFirst = gotSecond = NO;
		} 
	}
	if(ttt==0){
		i = [points count];
	NCITY = i;

	while (i--) {bPoint = (NXPoint *)[points elementAt:i];
  	   city[i].x=bPoint->x;city[i].y=bPoint->y;
	}
pathlength=0;
for ( i=0; i < NCITY; i++) {
     pathlength =pathlength+ get_distance( i, (i+1)%NCITY );
  } 
[pathL setIntValue:pathlength];
pathlength = pathlength*5;
  [initTemp setIntValue:pathlength];
  [freeze setIntValue:NCITY*10];}ttt=1;
	[self display];
	needsClearing = YES;
	return self;
}

int getSeparator(NXStream *stream)
/* 
 * Removes white space, commas, and plus signs from between 
 * number pairs.
 */
{
	int c, firstChar;
	
	NXUngetc(stream);
	c = firstChar = NXGetc(stream);
	while (NXIsSpace(c) || (c == ',') || (c == '+'))
		c = NXGetc(stream);
	/* 1 = ate something; 0 = ate nothing; -1 means EOF */
	if (c == firstChar)
		return 0;
	else if (c == EOF)
		return -1;
	NXUngetc(stream);
	return 1;
}

float getNumber(NXStream *stream)
/* 
 * Composes a floating point number from a string of number  
 * characters.
 */
{
	int c, i = 0;
	char temp[MAXNUMLENGTH];
	
	NXUngetc(stream);
	c = NXGetc(stream);
	while ((NXIsDigit(c) || (c == '.')) && (c != '\n') && (c != EOF)) {
		temp[i++] = c;
		c = NXGetc(stream);
	}
	NXUngetc(stream);
	temp[i] = 0;
	return (atof(temp));
}


static double get_distance( int k , int l )
{	
  static double distance;

  distance = sqrt( ( city[k].x - city[l].x ) * ( city[k].x - city[l].x )
			+ ( city[k].y - city[l].y ) * ( city[k].y - city[l].y ) );
			
  return distance;
}


static double get_uniform_random( )
{
  static double r;
	
  ir = (ialpha*ir) % ibeta;
	
  r = ( (double) ir ) / ( (double) ibeta );
  return r;
}


static double min_value ( double a, double b )
{
  if ( a < b ) return a;
  else return b;
}


@end