/*************************************************************************************************************************************/
/*                                                                                                                                   */
/*        C O N S T R U C T I O N   D ' U N E   C O U R B E   D E   V O N   K O C H  :                                               */
/*                                                                                                                                   */
/*                                                                                                                                   */
/*        Author of '$xtc/VonKoch.01$c' :                                                                                            */
/*                                                                                                                                   */
/*                    Jean-Francois Colonna (LACTAMME, AAAAMMJJhhmmss).                                                              */
/*                                                                                                                                   */
/*************************************************************************************************************************************/

#include  "INCLUDES.01.I"
                                        /* Introduit le 20051116095530...                                                            */

extern    void      *malloc();

extern    double    sqrt();

extern    double    cos();
extern    double    sin();

static    int       dimX=0;
#define   Xmin      0
#define   Xmax      (Xmin + (dimX-1))
                                        /* Definition des abscisses.                                                                 */
static    int       dimY=0;
#define   Ymin      0
#define   Ymax      (Ymin + (dimY-1))
                                        /* Definition des ordonnees.                                                                 */

#define   Xnormalise(x)                                                                                                                 \
                    ((int)((x-XG)/(XD-XG)*dimX))
#define   Ynormalise(y)                                                                                                                 \
                    ((int)((y-YB)/(YH-YB)*dimY))
#define   IMAGE(x,y)                                                                                                                    \
                    (*(image + ((int)(((y-Ymin)*dimX) + (x-Xmin)))))                                                                    \
                                        /* Acces a un point de l'image.                                                              */

#define   store(n,x,y)                                                                                                                  \
                    {                                                                                                                   \
                    if        (((x>=Xmin) && (x<=Xmax)) && ((y>=Ymin) && (y<=Ymax)))                                                    \
                              {                                                                                                         \
                              IMAGE(x,y) = n;                                                                                           \
                              }                                                                                                         \
                    else                                                                                                                \
                              {                                                                                                         \
                              }                                                                                                         \
                    }                                                                                                                   \
                                        /* Rangement d'un point valide d'une image.                                                  */

#define   XG        (0.00)
#define   XD        (1.00)
                                        /* Extremites Gauche et Droite sur l'axe des 'X'.                                            */
#define   YB        (0.00)
#define   YH        (1.00)
                                        /* Extremites Bas et Haut sur l'axe des 'Y'.                                                 */

#define   PRECIS    double                                                                                                              \
                                        /* Precision des calculs.                                                                    */

typedef   struct    complexe  {
                              PRECIS    reelle;
                              PRECIS    imaginaire;
                              }         complexe;
                                        /* 'struct' de definition d'un nombre complexe...                                            */
#define   Reelle(z)           (z.reelle)
#define   Imaginaire(z)       (z.imaginaire)

#define   Cinit(z,x1,y1)                                                                                                                \
                    {                                                                                                                   \
                    Reelle(z) = x1;                                                                                                     \
                    Imaginaire(z) = y1;                                                                                                 \
                    }
#define   Cegal(z,z1)                                                                                                                   \
                    {                                                                                                                   \
                    Cinit(z,Reelle(z1),Imaginaire(z1));                                                                                 \
                    }
#define   Csomme(z,z1,z2)                                                                                                               \
                    {                                                                                                                   \
                    complexe  zT;                                                                                                       \
                    Reelle(zT) = Reelle(z1) + Reelle(z2);                                                                               \
                    Imaginaire(zT) = Imaginaire(z1) + Imaginaire(z2);                                                                   \
                    Cegal(z,zT);                                                                                                        \
                    }
#define   Cdifference(z,z1,z2)                                                                                                          \
                    {                                                                                                                   \
                    complexe  zT;                                                                                                       \
                    Reelle(zT) = Reelle(z1) - Reelle(z2);                                                                               \
                    Imaginaire(zT) = Imaginaire(z1) - Imaginaire(z2);                                                                   \
                    Cegal(z,zT);                                                                                                        \
                    }
#define   Cproduit(z,z1,z2)                                                                                                             \
                    {                                                                                                                   \
                    complexe  zT;                                                                                                       \
                    Reelle(zT) = (Reelle(z1)*Reelle(z2)) - (Imaginaire(z1)*Imaginaire(z2));                                             \
                    Imaginaire(zT) = (Reelle(z2)*Imaginaire(z1)) + (Reelle(z1)*Imaginaire(z2));                                         \
                    Cegal(z,zT);                                                                                                        \
                    }
#define   Crotation(z,z1,rho,theta)                                                                                                     \
                    {                                                                                                                   \
                    complexe  zR;                                                                                                       \
                    Cinit(zR,rho*cos(theta),rho*sin(theta));                                                                            \
                    Cproduit(z,z1,zR);                                                                                                  \
                    }
#define   Cmodule2(z)                                                                                                                   \
                    ((Reelle(z)*Reelle(z)) + (Imaginaire(z)*Imaginaire(z)))

unsigned  char      *image;
                                        /* Definition de l'image a generer...                                                        */

PRECIS    seuil=0.002;
                                        /* Definition du seuil d'arret.                                                              */
PRECIS    angle=PI/3.0;
PRECIS    rapport=3.0;
                                        /* Definition de la transformation...                                                        */

#define   rotation_et_construction(theta)                                                                                               \
                    {                                                                                                                   \
                    Cegal(zGs3,zDs3);                                                                                                   \
                    Crotation(difference_reduite,difference_reduite,1.0,theta);                                                         \
                    Csomme(zDs3,zGs3,difference_reduite);                                                                               \
                    construction(zGs3,zDs3);                                                                                            \
                    }

construction(zG,zD)
complexe  zG,zD;
          {
          complexe  difference;

          Cdifference(difference,zD,zG);

          if        (sqrt(Cmodule2(difference)) < seuil)
                    {
                    store(BLANC,Xnormalise(Reelle(zG)),Ynormalise(Imaginaire(zG)));
                    store(BLANC,Xnormalise(Reelle(zD)),Ynormalise(Imaginaire(zD)));
                                        /* Marquage et arret de la recursivite...                                                    */
                    }
          else
                    {
                    complexe  zGs3,zDs3;
                    complexe  difference_reduite;

                    Cegal(zDs3,zG);
                    Crotation(difference_reduite,difference,1.0/rapport,0.0);

                    rotation_et_construction(+0.0*angle);
                    rotation_et_construction(+1.0*angle);
                    rotation_et_construction(+4.0*angle);
                    rotation_et_construction(+1.0*angle);
                    }
          return(OK);
          }

main()
          {
          complexe  zG,zD;
                                        /* Definition des extremites de la courbe.                                                   */
          int       x,y;
                                        /* Definition des coordonnees d'initialisation de l'image.                                   */

          Get(dimX,"dimX");
          Get(dimY,"dimY");
                                        /* Recuperation des dimensions en 'X' et en 'Y' de l'image a generer.                        */

          image=malloc(dimY*dimX);
                                        /* Definition de l'image a generer...                                                        */

          for       (y=Ymin ; y<=Ymax ; y++)
                    {
                    for       (x=Xmin ; x<=Xmax ; x++)
                              {
                              store(NOIR,x,y);
                                        /* Initialisation de l'image finale...                                                       */
                              }
                    }

          Cinit(zG,XG,(YB+YH)/2);
          Cinit(zD,XD,(YB+YH)/2);
                                        /* Initialisation des extremites de la courbe.                                               */

          construction(zG,zD);
                                        /* Construction de la courbe.                                                                */

          write(1,image,dimX*dimY);
                                        /* Sortie de l'image...                                                                      */
          }