/*******************************************************
* This program calculates the discreet primitive P(x)  *
* of a user defined real function F(x) in n points by  *
* using a classical Runge-Kutta method of order 4      *
* (version with graphic display of primitive P(x)).    *
* ---------------------------------------------------- *
* SAMPLE RUN:                                          *
* (Calculate 20 points of the primitive of F(x)=x      *
*  from x=-5 to x=5 with a finesse of 2. the value of  *
*  the primitive at x=-5 will be given as 25/2).       *
*                                                      *
*   Value of x begin: -5                               *
*   Value of x end  : 5                                *
*   Value of primitive in x begin: 12.5                *
*   Number of points to calculate: 50                  *
*   Number of intermediate points: 2                   *
*                                                      *
* (Exact values are yi = xi*xi/2).                     *
*                                                      *
* A graph of the primitive P(x) is displayed.          *
*                                                      *
*       C++ version in API style by J-P Moreau, Paris. *
*                      (www.jpmoreau.fr)               *
********************************************************
Used files: 
  gprimitiv.mak,
  basis.h (optionnal), graph2d.cpp, graph2d.h,
------------------------------------------------------*/
#include <stdio.h>
#include <math.h>
#include <basis.h>     // SQRT, REAL etc.   

#include <windows.h>

HDC  hdc;
RECT rect;
HPEN hpen, hpenOld;

//Functions used here of graph2D.cpp
void CourbeXY(int,int,float *,REAL,REAL);
void Legendes (int,char *,char *,char *,int,int,char *,int,int,char *);

//Function of fourier.cpp
void AnalyticFourierHn(REAL,REAL,int,REAL *,REAL *);

//"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;
  REAL dx,dy;
  if (ix2<ix1) {
	Swap(&ix1,&ix2); Swap(&iy1,&iy2);
  }
  il=(int) SQRT((ix2-ix1)*(ix2-ix1)+(iy2-iy1)*(iy2-iy1));
  dx=(REAL) (ix2-ix1)/il; 
  dy=(REAL) (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 application
REAL fa,xa,xb;
int    finesse,n;
float  Y[200];    //vector to store ordinates of P(x)

// given user defined function F(x)
REAL F(REAL x) {
  return x;
}

/******************************************************
* Calculate the discreet primitive of a real function *
* F(x) by using a Runge-Kutta method of order 4.      *
* --------------------------------------------------- *
* INPUTS:                                             *
*           a        begin x                          *
*           b        end x                            *
*           fa       value of primitive at x=a        *
*           n        number of calculated points      *
*           fi       number of intermediate points    *
*                    (finesse)                        *
*                                                     *
* OUTPUTS:                                            *
*           n ordinates put in vector Y(i).           *
* --------------------------------------------------- *
* Ref.: "Mathematiques en Turbo Pascal By Alain Rever-*
*        chon and Marc Ducamp, Armand Colin Editor,   *
*        Paris, 1991" [BIBLI 03].                     *
******************************************************/
void Primitive(REAL a,REAL b,REAL fa, int n, int fi) {
  REAL h,h2,h6,ly,lr,x,y,u;
  int i,j; long ni;
  h=(b-a)/fi/n; h2=h/2; h6=h/6;
  lr=fa; ly=F(a);
  Y[0]=(float) lr;
  for (i=1; i<=n; i++) {
    ni=(i-1)*fi-1;
    for (j=1; j<=fi; j++) {
      x=a+h*(ni+j);
      u=4*F(x+h2);
      x=x+h;
      y=F(x);
      lr=lr+h6*(ly+u+y); ly=y;
    }
    Y[i]=(float) lr; //current ordinate of primitive P(x)
  }
}

//calculate and draw discreet primitive of given function F(x)
void Exec_gprimitiv() {
  extern MaxX, MaxY;
  xa=-5; xb=5; fa=12.5; n=50; finesse=2;
  //calculate discreet primitive of F(x)
  Primitive(xa,xb,fa,n,finesse);
  //draw primitive P(x)
  MaxX=(rect.right-rect.left); //drawing zone in pixels
  MaxY=(rect.bottom-rect.top);
  CourbeXY(n+1,10,Y,xa,xb);
  Legendes(10,"  Primitive P(x)  ","X","Y",MaxX-220,50,"Curve: 50 points",
	                                       MaxX-220,80,"from x=-5 to x=5");
}

//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);
            hpen = CreatePen(PS_SOLID, 1, RGB(0, 0, 255));
            hpenOld = SelectObject(hdc,hpen);

            Exec_gprimitiv();
               
            SelectObject(hdc,hpenOld);
            DeleteObject(hpen);
            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[] = "Gprimitiv";
    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
        " Discreet primitive (analytic function)",    // 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 gprimitiv.cpp