/* Plot3DView.m   Copyright 1992 Steve Ludtke */
/* This view allows the display of 3d functions with real time rotation */

#import "Plot3DView.h"
#import "PControl.h"
#import <appkit/Control.h>
#import <appkit/Matrix.h>
#import <appkit/Application.h>
#import <appkit/NXColorWell.h>
#import <appkit/NXColorPanel.h>
#import <appkit/PrintInfo.h>
#import <dpsclient/event.h>
#import <dpsclient/psops.h>
#import <appkit/color.h>
#import <appkit/OpenPanel.h>
#import <stdio.h>
#import <math.h>
#include <libc.h>
#import "Expression.h"
#import "DensView.h"

void PSsethel();	/* makes Helvetica current font */

/* drawing and transformation macros */
#define MOVETO(x,y) pscom[psc]=dps_moveto; pspath[psc*2]=x; pspath[psc*2+1]=y; psc++ 
#define LINETO(x,y) pscom[psc]=dps_lineto; pspath[psc*2]=x; pspath[psc*2+1]=y; psc++
#define XFORM() x1=mx[0]*x+mx[1]*y+mx[2]*z; y1=mx[3]*x+mx[4]*y+mx[5]*z; z1=mx[6]*x+mx[7]*y+mx[8]*z

extern id NXApp;

@implementation Plot3DView
-initFrame:(NXRect *)myrect
{
int i;
[super initFrame:(NXRect *)myrect];
/* origin is in the center of the view, and size is +-1 + .1 for border */
[self setDrawSize:2.2 :2.2];
[self setDrawOrigin:-1.1 :-1.1];

/* bounding box for drawing */
psbbox[0]=psbbox[1]=-1.1;
psbbox[2]=psbbox[3]=1.1;

/* initialize preferences */
for (i=0; i<MAXSETS; i++) {
	pref[i].ndata=0; 
	pref[i].sym= -1;
	pref[i].color=NXConvertGrayToColor(0.0);
	pref[i].fileData=NULL;
	pref[i].expr=[[Expression alloc] init];
	[pref[i].expr parse:"0"];
}
/* initial function to display */
[pref[0].expr parse:"cos(sqrt(x^2+y^2))"];
pref[0].sym=6;

/* start spinning, all angles in radians */
chi=.4;
theta=.5236;
dchi=0.025;
Tgrids=20;
initflag=1;

[[[NXApp printInfo] setHorizCentered:YES] setOrientation:NX_LANDSCAPE andAdjust:YES];

timer=(DPSTimedEntry)DPSAddTimedEntry(TIMESTEP,itstime,self,NX_RUNMODALTHRESHOLD);

return self;
}

/* fix the scaling after being resized */
-superviewSizeChanged:(const NXSize *)oldsize
{
[super superviewSizeChanged:oldsize];
[self setDrawSize:2.2 :2.2];
[self setDrawOrigin:-1.1 :-1.1];
[self zoom:self];
return self;
}

-free
{
DPSRemoveTimedEntry(timer);
[super free];
return self;
}

/* please note that XFORM,MOVETO, and DRAWTO are macros, not functions */
/* when used, curly braces are necessary, since they are absent in the */
/* macro definition. */
-drawSelf:(NXRect *)myrect :(int)rectCount
{
static float mx[9];
char s[30];
float x,y,z,x1,y1,z1;
int i,j,k,fl,n;

n=Tgrids;

PSsetlinewidth(0.003);		/* not 0 since printer can't handle it */
				/* I'll come up with something better later */
				
/* calculate transformation matrix */
mx[0]=cos(chi);			/* theta=alt, chi=az  rotation matrix */
mx[1]= -cos(theta)*sin(chi);
mx[2]=sin(theta)*sin(chi);
mx[3]=sin(chi);
mx[4]=cos(theta)*cos(chi);
mx[5]= -sin(theta)*cos(chi);
mx[6]=0.0;
mx[7]=sin(theta);
mx[8]=cos(theta);

PSsetgray(NX_WHITE);		/* clear view */
NXRectFill(&bounds);

/* draw top and bottom of viewing area */
PSsetgray(.3333);
x=y=z=-0.5;
XFORM();
PSmoveto(x1,z1);
x=0.5;
XFORM();
PSlineto(x1,z1);
y=0.5;
XFORM();
PSlineto(x1,z1);
x=-0.5;
XFORM();
PSlineto(x1,z1);
y=-0.5;
XFORM();
PSlineto(x1,z1);
z=.5;
XFORM();
PSmoveto(x1,z1);
x=0.5;
XFORM();
PSlineto(x1,z1);
y=0.5;
XFORM();
PSlineto(x1,z1);
x=-0.5;
XFORM();
PSlineto(x1,z1);
y=-0.5;
XFORM();
PSlineto(x1,z1);
PSstroke();

psc=0;

/* plot the data */
for (i=0; i<MAXSETS; i++) {
	PSsetgray(NXGrayComponent(pref[i].color));
	if (pref[i].sym==-1) continue;
	switch(pref[i].sym) {
		case 0:		/* dot */
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			MOVETO(x1,z1);
			LINETO(x1,z1);
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
			}
		}
		break;
		case 1:		/* + */
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			MOVETO(x1-SS,z1);
			LINETO(x1+SS,z1);
			MOVETO(x1,z1-SS);
			LINETO(x1,z1+SS);
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
			}
		}
		break;
		case 2:		/* triangle */
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			MOVETO(x1-SS,z1-SS*.5773);
			LINETO(x1+SS,z1-SS*.5773);
			LINETO(x1,z1+SS*1.1547);
			LINETO(x1-SS,z1-SS*.5773);
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
			}
		}
		break;
		case 3:		/* box */
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			MOVETO(x1-SS,z1-SS);
			LINETO(x1+SS,z1-SS);
			LINETO(x1+SS,z1+SS);
			LINETO(x1-SS,z1+SS);
			LINETO(x1-SS,z1-SS);
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
			}
		}
		break;
		case 4:		/* X */
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			MOVETO(x1-SS,z1-SS);
			LINETO(x1+SS,z1+SS);
			MOVETO(x1-SS,z1+SS);
			LINETO(x1+SS,z1-SS);
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
			}
		}
		break;
		case 5:		/* 1/2 mesh */
		if (pref[i].fileData==NULL) fl=1;
		else fl=0;
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			if (j==0 || (x<pref[i].data[j-1].x&&pref[i].data[j-2].y==pref[i].data[j-1].y)) { MOVETO(x1,z1); }
			else { LINETO(x1,z1); }
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
				MOVETO(x1,z1);
			}
		}
		break;
		case 6:		/* mesh */
		if (pref[i].fileData!=NULL) break;
		for (j=0; j<pref[i].ndata; j++) {
			x=pref[i].data[j].x;
			y=pref[i].data[j].y;
			z=pref[i].data[j].z;
			XFORM();
			if (j%n==0) { MOVETO(x1,z1); }
			else { LINETO(x1,z1); }
			if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
				psc=0;
				MOVETO(x1,z1);
			}
		}
		for (j=0; j<n; j++) {
			for (k=j; k<pref[i].ndata; k+=n) {
				x=pref[i].data[k].x;
				y=pref[i].data[k].y;
				z=pref[i].data[k].z;
				XFORM();
				if (k==j) { MOVETO(x1,z1); }
				else { LINETO(x1,z1); }
					if (psc>950) {
				DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
					psc=0;
					MOVETO(x1,z1);
				}
			}
		}
		break;
	}
	if (psc!=0) {
		DPSDoUserPath(pspath,psc*2,dps_float,pscom,psc,psbbox,dps_ustroke);
		psc=0;
	}
}

/* display alt and az in upper right corner */
PSsetgray(NX_BLACK);
PSsethel();
sprintf(s,"Alt:%6.2f  Az:%6.2f",theta*57.295787,chi*57.295787);
PSmoveto(.6,1.05);
PSshow(s);
PSstroke();

return self;
}	

/* recalculates and redisplays data */
-zoom:sender
{
float x0,x1,xs,y0,y1,ys;
float z0,z1,x,y,zz;
int i,j=0,k,n;
id tmp;

/* calculate minima, maxima and step sizes */
x0=minX;
x1=maxX;
y0=minY;
y1=maxY;

Tgrids=n=[grids intValue];
if (n<=0||n>=60) { [grids setIntValue:20]; n=20; return self; }
xs=(x1-x0)/(float)n;
ys=(y1-y0)/(float)n;
x1-=xs/2.0;
y1-=ys/2.0;

/* make sure there's enough space in the data array */
for (i=0; i<NFN; i++) {
	if (pref[i].sym==-1 || pref[i].fileData!=NULL) continue;
	if (pref[i].ndata!=(n*n)) {
		if (pref[i].ndata!=0) free(pref[i].data);
		pref[i].data=malloc(sizeof(Point)*(n*n+2));
		pref[i].ndata=n*n;
	}
}
z0=MAXFLOAT;
z1= -MAXFLOAT;

/* calculate (function mode) or clip (file mode) the data */
/* also finds min/max Z */
for (i=0; i<NFN; i++) {
	if (pref[i].sym==-1) continue;
	tmp=pref[i].expr;
	if (pref[i].fileData!=NULL) {
		k=0;
		for (j=0; j<pref[i].nfdata; j++) {
			if (pref[i].fileData[j].x>x1||pref[i].fileData[j].x<x0||
			    pref[i].fileData[j].y>y1||pref[i].fileData[j].y<y0) continue;
			pref[i].data[k].x=(pref[i].fileData[j].x-x0)/(x1-x0)-.5;
			[tmp setVar:"x" value:pref[i].data[k].x];
			pref[i].data[k].y=(pref[i].fileData[j].y-y0)/(y1-y0)-.5;
			[tmp setVar:"y" value:pref[i].data[k].y];
			[tmp setVar:"z" value:pref[i].fileData[j].z];
			pref[i].data[k].z=[tmp resultValue];
			if (pref[i].data[k].z>z1) z1=pref[i].data[k].z;
			if (pref[i].data[k].z<z0) z0=pref[i].data[k].z;
			k++;
		}
		pref[i].ndata=k;
		continue;
	}
	
	j=0;
	for (y=y0; y<y1; y+=ys) {
		for (x=x0; x<x1; x+=xs) {
			pref[i].data[j].x=(x-x0)/(x1-x0)-.5;
			pref[i].data[j].y=(y-y0)/(y1-y0)-.5;
			[tmp setVar:"x" value:x];
			[tmp setVar:"y" value:y];
			zz=pref[i].data[j].z=[tmp resultValue];
			if (zz!=zz || zz<-MAXFLOAT || zz>MAXFLOAT) zz=pref[i].data[j].z=MAXFLOAT;
			if (zz>z1) z1=zz;
			if (zz<z0) z0=zz;
			j++;
		}
	}
}

if (z1<=z0) { z1=z0+.001; z0-=.001; }
/* allow the controller to override min/max Z */
[controller minmaxZ:&z0 :&z1];

/* scale and clip Z */
for (i=0; i<NFN; i++) {
	if (pref[i].sym==-1) continue;
	for (k=0; k<pref[i].ndata; k++) {
		if (pref[i].data[k].z>z1) pref[i].data[k].z=z1;
		if (pref[i].data[k].z<z0) pref[i].data[k].z=z0;
		pref[i].data[k].z=(pref[i].data[k].z-z0)/(z1-z0)-.5;
	}
}

/* display 3d plot */
[self display];
/* display density plot */
[controller updDen];
return self;
}

/* obvious */
-setAng:(float)az :(float)alt
{
theta=alt;
chi=az;
return self;
}

/* allows spinning and zooming with the mouse */
-mouseDown:(NXEvent *)oevent 
{
int oldMask,loop=1;
float ix=0,iy=0,ix1=0,iy1=0,ix2=0,iy2=0,ix3=0,iy3=0;
float mx[9],t;
NXEvent *event,evs;
long tm,tm2;

evs=*oevent;
oevent=&evs;
[self convertPoint:&oevent->location fromView:nil];
ix2=ix=oevent->location.x;
iy2=iy=oevent->location.y;
tm2=tm=oevent->time;

if (zMode) {
	mx[0]=cos(chi)/2.0;		/* theta=alt, chi=az  rotation matrix */
	mx[1]= -cos(theta)*sin(chi)/2.0;
	mx[2]=sin(theta)*sin(chi);
	mx[3]=sin(chi);
	mx[4]=cos(theta)*cos(chi);
	mx[5]= -sin(theta)*cos(chi);
	mx[6]=0.0;
	mx[7]=sin(theta)/2.0;
	mx[8]=cos(theta);
	ix=mx[0]*ix2+mx[1]*iy2+.5*mx[2];
	iy=mx[6]*ix2+mx[7]*iy2+.5*mx[8];
	ix1=ix2=ix3=ix;
	iy1=iy2=iy3=iy;
}

oldMask = [window addToEventMask:NX_LMOUSEDRAGGEDMASK];

while (loop) {
    event = [NXApp getNextEvent:(NX_LMOUSEUPMASK | NX_LMOUSEDRAGGEDMASK)];
    [self convertPoint:&event->location fromView:nil];

    if (zMode) {
        switch (event->type) {
        case NX_LMOUSEUP:
                loop = 0;
		ix=minX;
		iy=minY;
		ix1=maxX-ix;
		iy1=maxY-iy;
		ix2=oevent->location.x;
		iy2=oevent->location.y;
		ix3=event->location.x;
		iy3=event->location.y;
		if (ix3<ix2) { t=ix3; ix3=ix2; ix2=t; }
		if (iy3<iy2) { t=iy3; iy3=iy2; iy2=t; }
		ix2=ix2*ix1/2.0+ix1/2.0+ix;
		ix3=ix3*ix1/2.0+ix1/2.0+ix;
		iy2=iy2*iy1/2.0+iy1/2.0+iy;
		iy3=iy3*iy1/2.0+iy1/2.0+iy;
		
		if (ix2<ix3 && iy2<iy3) {
			minX=ix2;
			minY=iy2;
			maxX=ix3;
			maxY=iy3;
		}
		[controller setMM:minX :maxX :minY :maxY];
		[self zoom:self];
            break;
        case NX_LMOUSEDRAGGED:
		[self lockFocus];
		PSsetlinewidth(0.003);
		PSsetgray(1.0);
		PSmoveto(ix,iy);
		PSlineto(ix1,iy1);
		PSlineto(ix2,iy2);
		PSlineto(ix3,iy3);
		PSlineto(ix,iy);
		PSstroke();
		ix1=mx[0]*oevent->location.x+mx[1]*event->location.y+.5*mx[2];
		iy1=mx[6]*oevent->location.x+mx[7]*event->location.y+.5*mx[8];
		ix2=mx[0]*event->location.x+mx[1]*event->location.y+.5*mx[2];
		iy2=mx[6]*event->location.x+mx[7]*event->location.y+.5*mx[8];
		ix3=mx[0]*event->location.x+mx[1]*oevent->location.y+.5*mx[2];
		iy3=mx[6]*event->location.x+mx[7]*oevent->location.y+.5*mx[8];
		PSsetgray(0.0);
		PSmoveto(ix,iy);
		PSlineto(ix1,iy1);
		PSlineto(ix2,iy2);
		PSlineto(ix3,iy3);
		PSlineto(ix,iy);
		PSstroke();
		[self unlockFocus];
		[window flushWindow];
	    break;
	}
    }
    else {
        switch (event->type) {
        case NX_LMOUSEUP:
                loop = 0;
	        dchi=(event->location.x-ix2)/(float)(event->time-tm2)*10.0;
		if (fabs(dchi)<.004) dchi=0.0;
                [self display];
            break;
        case NX_LMOUSEDRAGGED:
	        theta+=(event->location.y-iy)*3.0;
	        chi+=(event->location.x-ix)*4.0;
	        if (theta>M_PI/2.0) theta=M_PI/2.0;
	        if (theta<-M_PI/2.0) theta=-M_PI/2.0;
	        while (chi>2.0*M_PI) chi-=2.0*M_PI;
	        while (chi<0.0) chi+=2.0*M_PI;
	        ix2=ix;
	        iy2=iy;
	        tm2=tm;
	        ix=event->location.x;
	        iy=event->location.y;
	        tm=event->time;
	        [self display];
            break;
        default:
            break;
        }
    }
}
[window setEventMask:oldMask];
return self;
}

/* function called by timer */
void itstime(DPSTimedEntry entry,double now,id call)
{
[call step];
return;
}

/* do one time step */
-step
{
if (initflag) { 
	/* first time, do initialization stuff */
	if (controller==nil) return self;
	[controller startup:pref];
	[self zoom:self];
	initflag=0; 
}
if (dchi==0.0) return self;
chi+=dchi;
while (chi>2.0*M_PI) chi-=2.0*M_PI;
while (chi<0.0) chi+=2.0*M_PI;
[self display];
return self;
}

/* set mouse mode */
-setMode:sender
{
zMode=[[sender selectedCell] tag];
return self;
}

/* new min/max values */
-zoomTo:(float)minx :(float)miny :(float)maxx :(float)maxy
{
minX=minx;
minY=miny;
maxX=maxx;
maxY=maxy;
return self;
}

-(int)acceptsFirstMouse { return (YES); }

-setcontroller:con 
{
controller=con;
return self;
}

/* allow user to pause spinning */
-togFreeze:sender
{
static float Tdchi;

if ([sender intValue]) {
	Tdchi=dchi;
	dchi=0.0;
	return self;
}
dchi=Tdchi;
return self;
}
@end
