/*************************************************************************************************************************************/
/* */
/* M A X I M I S A T I O N D E L A P L U S P E T I T E D I S T A N C E */
/* A L ' I N T E R I E U R D ' U N N U A G E D E P O I N T S T R I D I M E N S I O N N E L */
/* E N C O O R D O N N E E S S P H E R I Q U E S */
/* P A R U N E M E T H O D E D U T Y P E " R E C U I T S I M U L E " */
/* A V E C T E M P E R A T U R E : */
/* */
/* */
/* Author of '$xtc/recuit_si.44$c' : */
/* */
/* Jean-Francois COLONNA (LACTAMME, AAAAMMJJhhmmss). */
/* */
/*************************************************************************************************************************************/
#include "INCLUDES.01.I"
/* Introduit le 20051116102508... */
#define nPOINTS \
8 \
/* Nombre de points du nuage... */
#define nITERATIONS \
90000 \
/* Nombre d'iterations maximum... */
#define RHO__MINIMUM \
(MIN2(dimX,dimY)/2)
#define RHO__MAXIMUM \
RHO__MINIMUM
#define RHO__DELTA_INITIAL \
0.0
#define RHO__DELTA_FINAL \
0.0
/* Definition de la coordonnee spherique 'rho' initialise tel que rho=constante a priori... */
#define THETA_MINIMUM \
0
#define THETA_MAXIMUM \
PI
#define THETA_DELTA_INITIAL \
0.5
#define THETA_DELTA_FINAL \
0.1
/* Definition de la coordonnee spherique 'theta'. */
#define PHI__MINIMUM \
THETA_MINIMUM
#define PHI__MAXIMUM \
(2*THETA_MAXIMUM)
#define PHI__DELTA_INITIAL \
(2*THETA_DELTA_INITIAL)
#define PHI__DELTA_FINAL \
(2*THETA_DELTA_FINAL)
/* Definition de la coordonnee spherique 'phi'. */
#define PROBABILITE_DE_PERTURBER_UN_POINT \
0.1 \
/* Probabilite de perturber le point courant... */
#define PROBABILITE_D_ACCEPTER_UNE_MAUVAISE_SOLUTION \
0.05 \
/* Probabilite d'accepter une mauvaise solution... */
#define TEMPS \
(((double)(iterations-1))/((double)(nITERATIONS-1))) \
/* Simulation d'un temps dans [0,1]. */
extern double drand48();
#define PERTURBATION_ALEATOIRE(delta_initial,delta_final) \
((drand48() - 0.5)*((delta_initial*(1-temps))+(delta_final*temps))) \
/* Plus grande variation possible des coordonnees {x,y}, en rappelant que : */ \
/* */ \
/* drand48() E [0,1] */ \
/* */ \
/* d'ou la translation de 0.5 qui permet d'avoir un deplacement algebrique dans n'importe */ \
/* quelle direction et ce de facon isotrope... */
typedef struct point {
double x;
double y;
double z;
double rho;
double theta;
double phi;
} point;
/* Structure de definition d'un point... */
#define X(P) \
(P.x)
#define Y(P) \
(P.y)
#define Z(P) \
(P.z)
#define Rho_(P) \
(P.rho)
#define Theta(P) \
(P.theta)
#define Phi_(P) \
(P.phi)
/* Definition de l'acces aux coordonnees d'un point. */
extern double cos();
extern double sin();
#define ConversionX(rho,theta,phi) \
((rho)*cos(phi)*sin(theta))
#define ConversionY(rho,theta,phi) \
((rho)*sin(phi)*sin(theta))
#define ConversionZ(rho,theta,phi) \
((rho)*cos(theta))
/* Conversion spheriques-cartesiennes. */
extern double sqrt();
#define Carre(x) \
((x)*(x)) \
/* Definition du carre de 'x'. */
#define DistanceEuclidienne(A,B) \
sqrt(Carre(X(B)-X(A)) + Carre(Y(B)-Y(A)) + Carre(Z(B)-Z(A))) \
/* Definition de la distance entre deux points 'A' et 'B'. */
#define PERTURBATION(coordonnee,delta_initial,delta_final,minimum,maximum,periodiser) \
{ \
int perturber=VRAI; \
double coordonnee_perturbee; \
double temps=TEMPS; \
\
while (perturber == VRAI) \
{ \
coordonnee_perturbee = coordonnee + PERTURBATION_ALEATOIRE(delta_initial,delta_final); \
/* Perturbation de la coordonnee courante... */ \
\
if (periodiser == VRAI) \
{ \
if (coordonnee_perturbee < minimum) \
{ \
coordonnee_perturbee = coordonnee_perturbee - minimum + maximum; \
} \
else \
{ \
} \
\
if (coordonnee_perturbee > maximum) \
{ \
coordonnee_perturbee = coordonnee_perturbee - maximum + minimum; \
} \
else \
{ \
} \
} \
else \
{ \
} \
\
if ((coordonnee_perturbee < minimum) || (coordonnee_perturbee > maximum)) \
{ \
/* Les coordonnees hors-ecran sont rejetees... */ \
} \
else \
{ \
coordonnee = coordonnee_perturbee; \
perturber=FAUX; \
/* On s'arrete sur la premiere coordonnee perturbee situee dans l'ecran... */ \
} \
} \
}
#define DISTANCE_MINIMALE(distance_minimale) \
{ \
distance_minimale=+1e+300; \
\
for (i=0 ; i<nPOINTS ; i++) \
{ \
for (j=0 ; j<nPOINTS ; j++) \
{ \
if (i > j) \
{ \
double distance_courante=DistanceEuclidienne(Lpoints[i],Lpoints[j]); \
distance_minimale=MIN2(distance_minimale,distance_courante); \
/* Recherche de la distance minimale courante du nuage de points... */ \
} \
else \
{ \
} \
} \
} \
}
extern void *malloc();
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. */
static int dimZ=0;
#define Zmin 0
#define Zmax (Zmin + (dimZ-1))
/* Definition des profondeurs. */
#define IMAGE(x,y) \
(*(image + (((((int)y)-Ymin)*dimX) + (((int)x)-Xmin)))) \
/* Acces a un point de l'image. */
#define PROJECT 0.0
#define store(n,x,y,z) \
{ \
int xp=(x+(dimX/2))+((int)(PROJECT*((double)(z+(dimZ/2))))); \
int yp=(y+(dimY/2))+((int)(PROJECT*((double)(z+(dimZ/2))))); \
int niveau=n; \
/* "Pseudo-projection" du point {x,y,z}. */ \
\
if (((xp>=Xmin) && (xp<=Xmax)) && ((yp>=Ymin) && (yp<=Ymax))) \
{ \
IMAGE(xp,yp) = MAX2(NOIR+1,MIN2(BLANC,n)); \
} \
else \
{ \
} \
} \
/* Rangement d'un point valide d'une image. */
#define initialisation(x,y) \
{ \
if (((x>=Xmin) && (x<=Xmax)) && ((y>=Ymin) && (y<=Ymax))) \
{ \
IMAGE(x,y) = NOIR; \
} \
else \
{ \
} \
} \
/* Initialisation de l'image. */
main()
{
int generer_les_positions_finales=VRAI;
/* Pour controler la generation des positions finales uniquement (=1) ou des positions */
/* finales et intermediaires (=0). */
int sortir_l_image=VRAI;
/* Pour controler la sortie (=1) ou pas (=0) de l'image... */
int x,y;
/* Definition des coordonnees. */
unsigned char *image;
/* Definition de l'image a generer... */
int i,j,n;
int iterations;
/* Index divers... */
point Lpoints[nPOINTS];
point SLpoints[nPOINTS];
/* Liste des points et sa Sauvegarde. */
double distance_minimale_avant;
double distance_minimale_apres;
int distance_minimale_apres_valide=FAUX;
Get(dimX,"dimX");
Get(dimY,"dimY");
Get(dimZ,"dimZ");
dimZ = 4*dimZ;
/* Recuperation des dimensions en 'X', en 'Y' et en 'Z' de l'image a generer. On notera la */
/* multiplication par 4 de 'dimZ' car il est en general 4 fois plus petit que les autres... */
image=malloc(dimY*dimX);
/* Definition de l'image a generer... */
for (y=Ymin ; y<=Ymax ; y++)
{
for (x=Xmin ; x<=Xmax ; x++)
{
initialisation(x,y);
/* Initialisation de l'image a generer... */
}
}
for (n=0 ; n<nPOINTS ; n++)
{
double aveugle1=drand48();
double aveugle2=drand48();
double aveugle3=drand48();
double aveugle4=drand48();
double aveugle5=drand48();
double aveugle6=drand48();
/* Necessaire (c'est l'experience qui le montre...). */
Rho_(Lpoints[n]) = RHO__MINIMUM;
Theta(Lpoints[n]) = THETA_MINIMUM;
Phi_(Lpoints[n]) = PHI__MINIMUM;
/* Initialisation arbitraire des points. L'experience montre qu'une configuration */
/* initiale tres fortement symetrique est preferable a une situation aleatoire... */
X(Lpoints[n]) = ConversionX(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
Y(Lpoints[n]) = ConversionY(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
Z(Lpoints[n]) = ConversionZ(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
/* Et conversion... */
}
for (iterations=1 ; iterations<=nITERATIONS ; iterations++)
{
for (n=0 ; n<nPOINTS ; n++)
{
Rho_(SLpoints[n]) = Rho_(Lpoints[n]);
Theta(SLpoints[n]) = Theta(Lpoints[n]);
Phi_(SLpoints[n]) = Phi_(Lpoints[n]);
X(SLpoints[n]) = X(Lpoints[n]);
Y(SLpoints[n]) = Y(Lpoints[n]);
Z(SLpoints[n]) = Z(Lpoints[n]);
/* Sauvegarde de l'etat courant des points... */
}
if (distance_minimale_apres_valide == FAUX)
{
if (iterations == 1)
{
DISTANCE_MINIMALE(distance_minimale_avant);
/* Recherche de la distance minimale avant la perturbation. */
}
else
{
/* Dans ce cas 'distance_minimale_avant' est encore valide (calculee lors d'une iteration */
/* anterieure...). */
}
}
else
{
distance_minimale_avant = distance_minimale_apres;
}
if (iterations == 1)
{
fprintf(stderr,"distance minimale initiale=%f\n",distance_minimale_avant);
}
else
{
}
for (n=0 ; n<nPOINTS ; n++)
{
if (drand48() < PROBABILITE_DE_PERTURBER_UN_POINT)
{
double Rho__perturbe=Rho_(Lpoints[n]);
double Theta_perturbe=Theta(Lpoints[n]);
double Phi__perturbe=Phi_(Lpoints[n]);
PERTURBATION(Rho__perturbe,RHO__DELTA_INITIAL,RHO__DELTA_FINAL,RHO__MINIMUM,RHO__MAXIMUM,FAUX);
PERTURBATION(Theta_perturbe,THETA_DELTA_INITIAL,THETA_DELTA_FINAL,THETA_MINIMUM,THETA_MAXIMUM,FAUX);
PERTURBATION(Phi__perturbe,PHI__DELTA_INITIAL,PHI__DELTA_FINAL,PHI__MINIMUM,PHI__MAXIMUM,FAUX);
Rho_(Lpoints[n]) = Rho__perturbe;
Theta(Lpoints[n]) = Theta_perturbe;
Phi_(Lpoints[n]) = Phi__perturbe;
/* Perturbation aleatoire des points. */
X(Lpoints[n]) = ConversionX(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
Y(Lpoints[n]) = ConversionY(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
Z(Lpoints[n]) = ConversionZ(Rho_(Lpoints[n]),Theta(Lpoints[n]),Phi_(Lpoints[n]));
/* Et conversion... */
}
else
{
}
}
DISTANCE_MINIMALE(distance_minimale_apres);
/* Recherche de la distance minimale apres la perturbation. */
if ( (distance_minimale_apres < distance_minimale_avant)
|| (drand48() < PROBABILITE_D_ACCEPTER_UNE_MAUVAISE_SOLUTION)
)
/* On cherche a maximiser la distance minimale. On notera le test '<' (alors que la */
/* logique voudrait '<=') a cause du cas ou les points occuperaient initialement tous */
/* la meme position (d'ou : distance_minimale_avant=0) : dans ce cas, comme un certain */
/* nombre de points ne sont pas perturbes, on a donc distance_minimale_apres=0, d'ou avec */
/* avec le test '<=' le rejet de cette solution ; ainsi, on risque d'etre bloque */
/* indefiniment puisque l'on restaure alors l'etat precedent... */
{
distance_minimale_apres_valide = FAUX;
/* L'etat anterieur est inchange... */
for (n=0 ; n<nPOINTS ; n++)
{
Rho_(Lpoints[n]) = Rho_(SLpoints[n]);
Theta(Lpoints[n]) = Theta(SLpoints[n]);
Phi_(Lpoints[n]) = Phi_(SLpoints[n]);
X(Lpoints[n]) = X(SLpoints[n]);
Y(Lpoints[n]) = Y(SLpoints[n]);
Z(Lpoints[n]) = Z(SLpoints[n]);
/* Restauration de l'etat anterieur lorsqu'il y a degradation des performances... */
}
}
else
{
distance_minimale_apres_valide = VRAI;
/* On conserve un etat perturbe qui maximise la distance minimale... */
}
if (generer_les_positions_finales == 1)
{
}
else
{
for (n=0 ; n<nPOINTS ; n++)
{
store((BLANC*iterations)/nITERATIONS,X(Lpoints[n]),Y(Lpoints[n]),Z(Lpoints[n]));
/* Marquage des points... */
}
}
if (iterations == nITERATIONS)
{
fprintf(stderr,"distance minimale finale..=%f\n",distance_minimale_avant);
/* ATTENTION : c'est la distance minimale a l'avant-derniere iteration qui est editee car, */
/* en effet, il n'est pas sur que 'distance_minimale_apres' etait maximale... */
}
else
{
}
}
if (generer_les_positions_finales == 1)
{
for (n=0 ; n<nPOINTS ; n++)
{
store((BLANC*(Z(Lpoints[n])+(dimZ/2))/dimZ),X(Lpoints[n]),Y(Lpoints[n]),Z(Lpoints[n]));
/* Marquage des points... */
}
}
else
{
}
if (sortir_l_image == 1)
{
write(1,image,dimX*dimY);
/* Sortie de l'image... */
}
else
{
}
}