/******************************************************
*                  3D  SURFACES                       *
* --------------------------------------------------- *
* This program draws 3D surfaces defined by equations *
* Z = F(X,Y) by using unit graph_3d.pas.              *
* --------------------------------------------------- *
* SAMPLE RUN:                                         *
* ( Draw the surface defined by:                      *
*   Z = 8 sin(sqrt(X*X+Y*Y))/sqrt(X*X+Y*Y) )          *
*                                                     *
*  Input intervals [x1,x2] and [y1,y2]:               *
*    X1 X2 = -10 10                                   *
*    Y1 Y2 = -10 10                                   *
*  Number of lines: 100                               *
*  Number of pts per line : 100                       *
*  Real or uniform view (r/u): u                      *
*                                                     *
*  Choice of projection type:                         *
*  =========================                          *
*    1. real perspective                              *
*    2. ordinary parallel                             *
*    3. dimetric parallel                             *
*    4. isometric parallel                            *
*                                                     *
*  Your choice (1 to 4): 4                            *
* --------------------------------------------------- *
* Ref.: "Graphisme dans le plan et dans l'espace en   *
*        Turbo Pascal 4.0 By R. Dony, MASSON Paris,   *
*        1990" [BIBLI 12].                            *
*                                                     *
*      C++ version in API style by J-P Moreau, Paris  *
*                (with use of vmblock.cpp).           *
*                    (www.jpmoreau.fr)                *
******************************************************/
#include <windows.h>
#include <basis.h>
#include <vmblock.h>

#define  MACHEPS 1e-10
#define  perspective  1
#define  parallele    2

HDC  hdc;
RECT rect;

//Functions used here of module Graph_3D

//define drawing zone (ix1,ix2,iy1,iy2) in screen pixels
void Cloture(int,int,int,int);
//draw a frame around drawing zone (optional) 
void Bordure();

void Projette(double,double,double);
void InitProj(void);
void MoveXYZ(double,double,double);
void LineXYZ(double,double,double);


//"home made" graphic commands for hdc environment used by above functions
void DrawPixel(int ix,int iy) {
  //sorry, no other available command found
  Rectangle(hdc,rect.left+ix,rect.top+iy,
	            rect.left+ix+2,rect.top+iy+1);
}

void Swap(int *i1,int *i2) {
  int it;
  it=*i1; *i1=*i2; *i2=it;
}

void DrawLine(int ix1,int iy1,int ix2,int iy2) {
  int i,il,ix,iy;
  float dx,dy;
  if (ix2<ix1) {
	Swap(&ix1,&ix2); Swap(&iy1,&iy2);
  }
  il=(int) sqrt((ix2-ix1)*(ix2-ix1)+(iy2-iy1)*(iy2-iy1));
  dx=(float) (ix2-ix1)/il; 
  dy=(float) (iy2-iy1)/il;
  ix=ix1; iy=iy1;
  for (i=1;i<il+1; i++) {
	DrawPixel(ix,iy);
	ix=(int) (ix1+i*dx); 
	iy=(int) (iy1+i*dy);
  }
}

void OutText(int ix,int iy,char *s) {
  RECT rct;
  rct.left=ix; rct.right=ix+10*strlen(s);
  rct.top=iy-5; rct.bottom=iy+15;
  DrawText(hdc, s, -1, &rct,
                DT_SINGLELINE | DT_CENTER | DT_VCENTER);
}

//*******************************************
//global variables of Surfaces program
int    *Hmax,*Hmin;
double x1,x2,yy1,y2;  //y1 reserved by system
double incx,incy;
double f1,f2,f3,f4;
int    c1,c2,c3,c4;
double echx, echy;
char   view;
int    xg,yg,xd,yd,nbrelignes,nbrepoints;

void   *vmblock;
int    LimX, LimY;

extern double de,phi,rho,theta;

//Surface Z=F(X,Y)
double F(double x, double y) {
  double k;
  k=sqrt(x*x+y*y);
  if (fabs(k)>MACHEPS)
    return (8*sin(k)/k);
  else
    return 8.0;
}

/*double F(double x, double y) {
  return (x*x+y*y);
} */

int Signe(double x) {
  if (x > 0.0) return 1;
  else if (x < 0.0) return -1;
  else return 0;
}

//data here are static (precompiled)
void Data() {
  extern projection;
  int choix;
  //define client drawing zone in screen pixels
  LimX=(rect.right-rect.left);
  LimY=(rect.bottom-rect.top);
  //choose range for x and y
  x1=-16; x2=16; yy1=-12; y2=12;
  //choose number of lines
  nbrelignes=80;
  //choose number of points per line
  nbrepoints=120;
  //choose type of view (real or uniform)
  view='u';
  //choose projection type
  // 1. real perspective
  // 2. ordinary parallel
  // 3. dimetric parallel
  // 4. isometric parallel
  choix=4;
  de=1000;
  projection=parallele;
  switch(choix) {
    case 1: rho=100; theta=30; phi=30; 
	        projection=perspective;
			break;
    case 2: theta=30; phi=30;
	        break;
    case 3: theta=22.20765; phi=20.704811;
            break;
    case 4: theta=45; phi=35.26439;
  }
}//Data()

void Init() {
  double aux; int i;
  incx=(x2-x1)/nbrepoints;
  incy=(y2-yy1)/nbrelignes;
  c1=1; c2=LimX; c3=1; c4=LimY;
  f1=1e10; f2=-f1; f3=f1; f4=-f1;
  xg=-1; yg=-1; xd=-1; yd=-1;
  for (i=0; i<=LimX; i++)  Hmax[i]=0;
  for (i=0; i<=LimX; i++)  Hmin[i]=LimY;
  if (theta < 0 || theta > 180) {
    aux=x1; x1=x2; x2=aux; incx=-incx;
    aux=yy1; yy1=y2; y2=aux; incy=-incy;
  }
} //Init()

void Fenetre() {
  double x,y,z;
  int ligne,point;
  extern double xproj,yproj; //declared in module graph3d.cpp
  for (ligne=0; ligne<=nbrelignes; ligne++) {
    y=y2-ligne*incy;
    for (point=0; point<=nbrepoints; point++) {
      x=x1+point*incx;
      z=F(x,y);
      Projette(x,y,z);
      if (xproj < f1)  f1=xproj;
      if (xproj > f2)  f2=xproj;
      if (yproj < f3)  f3=yproj;
      if (yproj > f4)  f4=yproj;
    }
  }
} //Fenetre()

void Echelles() {
  echx=(c2-c1)/(f2-f1);
  echy=(c4-c3)/(f4-f3);
  if (view=='r') 
    if (echx < echy)  echy=echx; else echx=echy;
}

int Max(int x1,int x2) {
  if (x1>x2) return x1; else return x2;
}

int Min(int x1,int x2) {
  if (x1<x2) return x1; else return x2;
}

int Check(int x1,int y1,int x2,int y2) {
  if (x1>0 && x1<LimX && x2>0 && x2<LimX && y1>0 && y1<LimY && y2>0 && y2<LimY)
     return 1;
  else
     return 0;
}

void Horizon(int x1,int y1,int x2,int y2) {
  int x,y,dx;
  double pente;
  dx=Signe(x2-x1);
  if (dx==0 || Check(x1,y1,x2,y2)==0)
    if (y2>=0 && y2<=LimY) {
      Hmax[x2+1]=Max(Hmax[x2],y2);
      Hmin[x2+1]=Min(Hmin[x2],y2);
    }
    else {
      Hmax[x2+1]=Hmax[x2];
      Hmin[x2+1]=Hmin[x2];
    }
  else {
    pente=(y2-y1)/(x2-x1);
    for (x=x2+1; x<=x1; x++) {
      y=(int) floor(pente*(x-x1)+y1);
      Hmax[x]=Max(Hmax[x],y);
      Hmin[x]=Min(Hmin[x],y);
    }
  }
} //Horizon()

void Visibilite(int x,int y, int *visi) {
  if (y < Hmax[x] && y > Hmin[x])
    *visi=0;
  else if 
    (y >= Hmax[x]) *visi=1; 
  else 
    *visi=-1;
}

void Inter1(int x1,int y1,int x2,int y2, int *Tabaux, int *xi, int *yi) {
  double  ct1,ct2,p1,p2,xii,yii;
  if (x2-x1==0 || Check(x1,y1,x2,y2)==0)  {
    xii=1.0*x2;
    yii=1.0*Tabaux[x2];
  }
  else {
    p1=(double) (y2-y1)/(x2-x1);
    p2=(double) (Tabaux[x2]-Tabaux[x1])/(x2-x1);
    if (fabs(p1) > 1e-10 && fabs(p1-p2) > 1e-10) {
      ct1=(double) (y1-p1*x1); 
      ct2=(double) (Tabaux[x1]-p2*x1);
      yii=(p1*ct2-p2*ct1)/(p1-p2);
      xii=(yii-ct1)/p1;
    }
    else {
      xii=1.0*x2;
      yii=1.0*y2;
    }
  }
  *xi=(int) floor(xii);
  *yi=(int) floor(yii);
} //inter1

void AreteFermeture(int x,int y,int *xlateral,int *ylateral) {
  int ix,iy;
  ix=*xlateral; iy=*ylateral;
  if (*xlateral != -1)  Horizon(ix,iy,x,y);
  *xlateral=x;
  *ylateral=y;
}

//Hidden lines not removed
void DessinFonct() {
  double x, y, z;
  int ligne, point;
  int xcour,ycour;
  extern double xproj,yproj;
  extern xscreen, yscreen;
  for (ligne = 0; ligne <=nbrelignes; ligne++) {
    y = y2 - ligne * incy;
    x = x1;
    z = F(x, y);
    Projette(x, y, z);
    xcour = (int) (floor((xproj - f1) * echx) + c1);
    ycour = (int) (floor((yproj - f3) * echy) + c3);
    for (point = 0; point<=nbrepoints; point++) {
      x = x1 + point * incx;
      z = F(x, y);
      Projette(x, y, z);
      xscreen = (int) (floor((xproj - f1) * echx) + c1);
      yscreen = (int) (floor((yproj - f3) * echy) + c3);
      DrawLine(xcour, LimY - ycour, xscreen, LimY - yscreen);
      xcour = xscreen; ycour = yscreen;
    }
  }
} //DessinFonct()

//Hidden lines removed
void DessinFonction() {
  int ligne,point,xi,yi;
  int xprec,yprec,xcour,ycour;
  int visicour,visiprec;
  double x,y,z;
  extern double xproj,yproj;
  for (ligne=0; ligne<=nbrelignes; ligne++) {
    y=y2-ligne*incy;
    x=x1;
    z=F(x,y);
    Projette(x,y,z);
    xprec=(int) (floor((xproj-f1)*echx)+c1);
    yprec=(int) (floor((yproj-f3)*echy)+c3);
    AreteFermeture(xprec,yprec,&xd,&yd);
    Visibilite(xprec,yprec,&visiprec);
    for (point=0; point<=nbrepoints; point++) {
      x=x1+point*incx;
      z=F(x,y);
      Projette(x,y,z);
      xcour=(int) floor((xproj-f1)*echx)+c1;
      ycour=(int) floor((yproj-f3)*echy)+c3;
      Visibilite(xcour,ycour,&visicour);
      if (Hmax[xcour]==0 || Hmin[xcour]==LimY)
        visicour=visiprec;
      if (visicour==visiprec) {
        if (visicour==1 || visicour==-1) {
	      if (Check(xprec,LimY-yprec,xcour,LimY-ycour)) {
            DrawLine(xprec,LimY-yprec,xcour,LimY-ycour);
            Horizon(xprec,yprec,xcour,ycour);
          }
        }
      }
      else {
        if (visicour==0) {
          if (visiprec==1)
            Inter1(xprec,yprec,xcour,ycour,Hmax,&xi,&yi);
          else
            Inter1(xprec,yprec,xcour,ycour,Hmin,&xi,&yi);
          if (Check(xprec,LimY-yprec,xi,LimY-yi)) { 
            DrawLine(xprec,LimY-yprec,xi,LimY-yi);
            Horizon(xprec,yprec,xi,yi);
          }
        }
        else {
          if (visicour==1) {
            if (visiprec==0) {
              Inter1(xprec,yprec,xcour,ycour,Hmax,&xi,&yi);
              if (Check(xi,LimY-yi,xcour,LimY-ycour)) {
                DrawLine(xi,LimY-yi,xcour,LimY-ycour);
                Horizon(xi,yi,xcour,ycour);
              }
            }
            else {
              Inter1(xprec,yprec,xcour,ycour,Hmin,&xi,&yi);
	          if (Check(xprec,LimY-yprec,xi,LimY-yi)) {
                DrawLine(xprec,LimY-yprec,xi,LimY-yi);
                Horizon(xprec,yprec,xi,yi);
              }
              Inter1(xprec,yprec,xcour,ycour,Hmax,&xi,&yi);
			  if (Check(xi,LimY-yi,xcour,LimY-ycour)) {
                DrawLine(xi,LimY-yi,xcour,LimY-ycour);
                Horizon(xi,yi,xcour,ycour);
              }
            }
          }
          else {
            if (visiprec==0) {
              Inter1(xprec,yprec,xcour,ycour,Hmin,&xi,&yi);
              if (Check(xi,LimY-yi,xcour,LimY-ycour)) {
                DrawLine(xi,LimY-yi,xcour,LimY-ycour);
                Horizon(xi,yi,xcour,ycour);
              }
            }
            else {
              Inter1(xprec,yprec,xcour,ycour,Hmax,&xi,&yi);
	          if (Check(xprec,LimY-yprec,xi,LimY-yi)) {
                DrawLine(xprec,LimY-yprec,xi,LimY-yi);
                Horizon(xprec,yprec,xi,yi);
              }
              Inter1(xprec,yprec,xcour,ycour,Hmin,&xi,&yi);
	          if (Check(xi,LimY-yi,xcour,LimY-ycour)) {
                DrawLine(xi,LimY-yi,xcour,LimY-ycour);
                Horizon(xi,yi,xcour,ycour);
              }
            }
          }
        }
      } //if visicour==0 (else)
      visiprec=visicour;
      xprec=xcour;
      yprec=ycour;
    } //point loop
    AreteFermeture(xcour,ycour,&xg,&yg);
  } //ligne loop
} //DessinFonction()

void Affiche() {
  char aux[20],ch[100];
  extern projection;
  sprintf(aux,"X=[%4.1f", x1); strcpy(ch,aux);
  sprintf(aux,",%4.1f]",x2); strcat(ch,aux);
  sprintf(aux," Y=[%4.1f", yy1); strcat(ch,aux);
  sprintf(aux,",%4.1f]",y2); strcat(ch,aux);
  if (projection==perspective) {
    sprintf(aux," Rho=%4.1f",rho); strcat(ch,aux);
  }
  sprintf(aux," Theta=%5.3f",theta); strcat(ch,aux);
  sprintf(aux," Phi=%5.3f",phi); strcat(ch,aux);
  OutText(40,5,ch);
  sprintf(aux,"(%d lines,",nbrelignes); strcpy(ch,aux);
  sprintf(aux,"%d pts per line)",nbrepoints); strcat(ch,aux);
  OutText(155,20,ch);
  sprintf(ch,"f1=%10.2f f2=%10.2f f3=%10.2f f4=%10.2f",f1,f2,f3,f4);
  OutText(40,35,ch);
}

void Surfaces3D()  {
  Data();
  //Allocate memory space for Hmin,Hmax and Tabaux
  vmblock = vminit();
  Hmin   = (int *) vmalloc(vmblock, VEKTOR, LimX+1, 0);
  Hmax   = (int *) vmalloc(vmblock, VEKTOR, LimX+1, 0);
  if (! vmcomplete(vmblock)) {
    LogError ("Memory full!", 0, __FILE__, __LINE__);
    return;
  }	      
  Init();
  InitProj();
  Fenetre();
  Echelles();
  DessinFonction();
  vmfree(vmblock); //free memory
  Affiche();
}

//Handle window Paint and Destroy messages
long FAR PASCAL /*_export*/ WndProc(HWND hwnd,
    UINT message, UINT wParam, LONG lParam) 
{
    PAINTSTRUCT ps;
    
    switch (message) {
        case WM_PAINT:
			hdc = BeginPaint(hwnd, &ps);
			GetClientRect(hwnd, &rect);

            Surfaces3D();
               
            EndPaint(hwnd, &ps);
            return 0;
        case WM_DESTROY:
            PostQuitMessage(0);
            return 0;
    }
    return DefWindowProc(hwnd, message, wParam, lParam);
}

// main program
int PASCAL WinMain(HANDLE hInstance, HANDLE hPrevInstance,
    LPSTR lpszCmdParam, int nCmdShow)
{
    static char szAppName[] = "Surfaces 3D";
    HWND        hwnd;
    MSG         msg;
    WNDCLASS    wndclass;

    if (!hPrevInstance) {
        wndclass.style         = CS_HREDRAW | CS_VREDRAW;
        wndclass.lpfnWndProc   = WndProc;
        wndclass.cbClsExtra    = 0;
        wndclass.cbWndExtra    = 0;
        wndclass.hInstance     = hInstance;
        wndclass.hIcon         = LoadIcon(NULL, IDI_APPLICATION);
        wndclass.hCursor       = LoadCursor(NULL, IDC_ARROW);
        wndclass.hbrBackground = GetStockObject(WHITE_BRUSH);
        wndclass.lpszMenuName  = NULL;
        wndclass.lpszClassName = szAppName;

        RegisterClass(&wndclass);
    }

    hwnd = CreateWindow(szAppName,  // window class name
        " SURFACES 3D",             // window caption
        WS_OVERLAPPEDWINDOW,        // window style                       
        CW_USEDEFAULT,              // initial x position
        CW_USEDEFAULT,              // initial y position
        CW_USEDEFAULT,              // initial x size
        CW_USEDEFAULT,              // initial y size
        NULL,                       // parent window handle
        NULL,                       // window menu handle
        hInstance,                  // program instance handle
        NULL);                      // creation parameters

   ShowWindow(hwnd, nCmdShow);
   UpdateWindow(hwnd);

   while (GetMessage(&msg, NULL, 0, 0)) {
       TranslateMessage(&msg);
       DispatchMessage(&msg);
   }
   return msg.wParam;
}

//end of file surfaces.cpp