/**************************************************************************/
/* analog apfel3D.c nur mit Farbinterpretation der Flaechen               */
/**************************************************************************/
/* Apfel3d2.c                                                  28.12.87   */
/* Version 1.0                                                            */
/*                                                                        */
/*                                                                        */
/*                                                                        */
/* fuer Markt+Technik geschrieben von Heiko Knappe                        */
/*                                                                        */
/* Komplexe Dynamische Systeme                                            */
/**************************************************************************/


#include <exec/types.h>
#include <math.h>
#include <intuition/intuition.h>
#include <graphics/gfxmacros.h>

#define WIDTH   320  /* Screen Modus */
#define HEIGHT  250
#define DEPTH   5
#define FLAGS   0
#define FWIDTH  200
#define FHEIGHT 100
#define FXOFF   110
#define FYOFF   150

#define N 40         /* Anzahl der Stuetzpunkte in x,y Richtung */

/* falls das Bild gespeichert werden soll! */
/* extern short save(); */

struct IntuitionBase *IntuitionBase;
struct GfxBase *GfxBase;
struct Screen *MyScreen;
struct RastPort *MyRastPort;
struct AreaInfo MyAreaInfo;
struct TmpRas MyTmpRas;
struct NewScreen MyNewScreen=
 {0,0,WIDTH,HEIGHT,DEPTH,1,0,FLAGS,CUSTOMSCREEN,0,0,0,0};

int   n,nmax=100;            /* Iterationsschritte                  */
float emin=.001,emax=100000; /* Abbruchkriterium                    */
float rmin=-1.5,rmax=2.5;    /* Grenzen des Bildes                  */
float imin=-1.2,imax=1.2;
float rst,ist;
float re,im,re1,re2,im1,im2,cre,cim;
                             /* re,im laufender Re- und Im-Teil     */
                             /* re1,re2,im1,im2  -> Folgenglieder   */
float d;                     /* Abstand zu dem Bezugspunkt          */
int   i,j,k;
short land[N][N];            /* Feld mit den Punktwerten            */
short x[4],y[4];             /* 4 Eckpunkte der Flaechen            */
short stx,sty;
short anzcol,col1,col2,col;  /* Farbanzahl, Farbregister            */
short bx=1,by=1,bz=1;        /* Beleuchtungsvektor                  */
float b;
UWORD buffer[5*10];          /* Variablen fuer AreaFill             */
long  *speicher;


/***********************/
/* Hauptprogramm       */
/***********************/

main(argc,argv)
int argc;
char *argv[];
{
 if((argc!=2)&&(argc!=1))
 {
    printf("Eingabe: %s <nmax,rmin,rmax,imin,imax>\n",argv[0]);
    printf(" nmax: maximale Anzahl der Duerchlaeufe\n");
    printf(" [rmin,rmax] [imin,imax]: Wertebereich der Bildebene\n");
    close_things();
 }
 if(argc==2)
    if(sscanf(argv[1],"%d,%f,%f,%f,%f",&nmax,&rmin,&rmax,&imin,&imax)!=5)
    {
       printf("falsche Eingabe\n");
       close_things();
    }
 printf("\nParameter:\n emin=%f emax=%f\n",emin,emax);
 printf(" nmax=%d\n [%f,%f] [%f,%f]\n",nmax,rmin,rmax,imin,imax);
 emin*=emin;  /* mit dem Quadrat zu vergleichen erspart Wurzeln zu */
 emax*=emax;  /* berechnen -> Schnelligkeit                        */
 b=sqrt((float)(bx*bx+by*by+bz*bz));
 apfel();
 open_things();
 ausgabe();
 while(MyScreen->TopEdge<100);
 /* falls das Bild gespeichert werden soll      */
 /* save("Apfel.pic",1,MyScreen); */
 close_things();
}


/***********************/
/* Initialisierung     */
/***********************/

open_things()
{
 if((IntuitionBase=OpenLibrary("intuition.library",0))==0)
    close_things();
 if((GfxBase=OpenLibrary("graphics.library",0))==0)
    close_things();
 if((MyScreen=OpenScreen(&MyNewScreen))==0)
    close_things();
 MyRastPort=&(MyScreen->RastPort);
 init_colors();
 init_fill();
} 

init_fill()
{
 if((speicher=AllocRaster(WIDTH,HEIGHT))==0)
    close_things();
 MyRastPort->AreaInfo=&MyAreaInfo;
 MyRastPort->TmpRas=InitTmpRas(&MyTmpRas,speicher,RASSIZE(WIDTH,HEIGHT));
 InitArea(&MyAreaInfo,&buffer[0],10);
}

close_things()
{
 if(speicher!=0) FreeRaster(speicher,WIDTH,HEIGHT);
 if(IntuitionBase!=0) CloseLibrary(IntuitionBase);
 if(GfxBase!=0) CloseLibrary(GfxBase);
 if(MyScreen!=0) CloseScreen(MyScreen);
 exit();
} 

init_colors()
{
 for(i=0,anzcol=1;i<DEPTH;anzcol*=2,++i);
 SetRGB4(&(MyScreen->ViewPort),0,0,0,0); /* Hintergrund schwarz */
 for(i=1;i<anzcol/2;++i)
    SetRGB4(&(MyScreen->ViewPort),i,
            (short)((float)i/anzcol*15*2),1,0);
 for(i=anzcol/2;i<anzcol;++i)
    SetRGB4(&(MyScreen->ViewPort),i,15,1,
            (short)((float)(i-anzcol/2)/anzcol*15*2));
}


/***********************/
/* Berechnung der      */
/* Stuetzpunkte        */
/***********************/

apfel()
{
 printf("Berechnung der Stuetzpunkte \n");
 rst=(rmax-rmin)/N;
 ist=(imax-imin)/N;
 for(i=0;i<N;++i) /* in x-Richtung */
 {
    re=rmin+i*rst;
    for(j=0;j<N;++j) /* in y-Richtung */
    {
       im=imin+j*ist;
       for(n=re1=im1=0,cre=re,cim=im;n<nmax;++n)
       {
          im2=2*re1*im1-cim,re2=re1*re1-im1*im1-cre;
          im1=im2,re1=re2;
          d=im1*im1+re1*re1;
          /* Ist der Attraktor erreicht?                       */
          if(d<emin)  break;
          /* Attraktor Unendlich?                              */
          if(d>emax)  break;
       }
       land[i][j]=n;
    }
 }
 printf("-> fertig\n");
}


/***********************/
/* Ausgabe             */
/***********************/

ausgabe()
{
 stx=FWIDTH/N,sty=FHEIGHT/N;
 for(i=0;i<N-1;++i)
 {
     for(j=0;j<N-1;++j)
     {
        /* Berechnung der Eckpunkte der Flaechen */
        x[0]=FXOFF+stx*i-sty*j,y[0]=FYOFF-land[i][j]+sty*j;
        x[1]=FXOFF+stx*(i+1)-sty*j,y[1]=FYOFF-land[i+1][j]+sty*j;
        x[2]=FXOFF+stx*(i+1)-sty*(j+1),y[2]=FYOFF-land[i+1][j+1]+sty*(j+1);
        x[3]=FXOFF+stx*i-sty*(j+1),y[3]=FYOFF-land[i][j+1]+sty*(j+1);
        for(k=0;k<3;++k)
           if(y[k]<0) y[k]=0;
        zeich();
     }
 }
}

zeich()
{
 /* im Vergleich zu apfel3d2.c wird das entstandene Viereck noch in 2   */
 /* Dreiecke unterteilt um so richtige Farbberechnung zu gewaehrleisten */
 farbe();
 col=(col1+col2)/2;
 SetAPen(MyRastPort,anzcol-col);  /* Flaeche1 zeichnen */
 SetOPen(MyRastPort,col);
 AreaMove(MyRastPort,x[0],y[0]); AreaDraw(MyRastPort,x[1],y[1]);
 AreaDraw(MyRastPort,x[2],y[2]); AreaDraw(MyRastPort,x[3],y[3]);
 AreaEnd(MyRastPort);
}

farbe()  /* Berechnung der Farbe einer Flaeche */
{
 float n,h;
 int nx,ny,nz;
 short v1x,v1y,v1z,v2x,v2y,v2z;

 /* 1. Flaeche */
 v1x=1,v1y=0,v1z=land[i+1][j]-land[i][j];
 v2x=0,v2y=1,v2z=land[i+1][j+1]-land[i][j];
 nx=v1y*v2z-v1z*v2y,ny=v1z*v2x-v1x*v2z,nz=v1x*v2y-v1y*v2x;
 n=sqrt((float)(nx*nx+ny*ny+nz*nz));
 if(n==0) col1=1;
 h=acos((float)(nx*bx+ny*by+nz*bz)/b/n)-PI/2;
 h=(h<0)? -h : h;
 if((col1=(short)((PI/2-h)*(anzcol-1)*2/PI)+1)>(anzcol-1)) col1=anzcol-1;
 if(col1<=0) col1=1;
 /* 2. Flaeche */
 v1x=1,v1y=0,v1z=land[i][j+1]-land[i][j];
 v2x=0,v2y=1,v2z=land[i+1][j+1]-land[i][j];
 nx=v1y*v2z-v1z*v2y,ny=v1z*v2x-v1x*v2z,nz=v1x*v2y-v1y*v2x;
 n=sqrt((float)(nx*nx+ny*ny+nz*nz));
 if(n==0) col2=1;
 h=acos((float)(nx*bx+ny*by+nz*bz)/b/n)-PI/2;
 h=(h<0)? -h : h;
 if((col2=(short)((PI/2-h)*(anzcol-1)*2/PI)+1)>(anzcol-1)) col2=anzcol-1;
 if(col2<=0) col1=1;
}
