/******************************************************************
 *                                                                *
 *         Nest version 2 - 3D Wasp Nest Simulator                *
 *            Copyright (C) 1997 Sylvain GUERIN                   *
 *        LIASC, ENST Bretagne & Santa Fe Institute               *
 *                                                                *
 ******************************************************************
 
 E-mail: Sylvain.Guerin@enst-bretagne.fr
 or: Sylvain Guerin
 13 rue des Monts Lorgeaux, 51460 L'EPINE, FRANCE
 
 ******************************************************************
 
 This program is free software; you can redistribute it and/or
 modify it under the terms of the GNU General Public License
 as published by the Free Software Foundation; either version 2
 of the License, or (at your option) any later version.
 
 This program is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 GNU General Public License for more details.
 
 You should have received a copy of the GNU General Public License
 along with this program; if not, write to the Free Software
 Foundation, Inc., 59 Temple Place - Suite 330, 
 Boston, MA  02111-1307, USA.
 
 ******************************************************************
     
      Name          : pattern.c
      Component     : fitness evaluation
      Fonctionality : contains methods for class Pattern
  
 ******************************************************************/

#include "graph.h"
#include "pattern.h"
#include <math.h>

/* NOTA: pour le moment, dans les patterns, RED=BLUE */
/* Sylvain GUERIN, 15 mars 1998 */


/*****************/
/* Objet Pattern */
/*****************/

Pattern::Pattern (int aMaxSize, int anInitRuleId, SimulationType aSimuType)
{
  int i,j,k;

  simuType = aSimuType;
  initRuleId = anInitRuleId;
  maxSize = aMaxSize;
  cell = new Cell** [aMaxSize];
  for (i=0; i<maxSize; i++) {
    cell[i] = new Cell* [aMaxSize];
  }
  for (i=0; i<maxSize; i++) {
    for (j=0; j<maxSize; j++) {
      cell[i][j] = new Cell [aMaxSize];
    }
  }
  for (i=0; i<maxSize; i++) {
    for (j=0; j<maxSize; j++) {
      for (k=0; k<maxSize; k++) {
	cell[i][j][k].set (EMPTY, -1, -1);
      }
    }
  }
  
  minX = maxSize;
  minY = maxSize;
  minZ = maxSize;
  maxX = -1;
  maxY = -1;
  maxZ = -1;
  numberOfCells = 0;
  isEncoded = 0;
  classIndex = -1;
}

Pattern::~Pattern ()
{
  int i, j;
  
  for (i=0; i<maxSize; i++) {
    for (j=0; j<maxSize; j++) {
      delete[] cell[i][j];
    }
  }
  for (i=0; i<maxSize; i++) {
    delete[] cell[i];
  }
  delete[] cell;

  if (isEncoded) {
    freeEncoding ();
  }
}

void Pattern::setClassIndex (int anIndex)
{
  classIndex = anIndex;
}

int Pattern::getClassIndex ()
{
  return classIndex;
}

void Pattern::set (int x, int y, int z, 
		   typeOfCell aCellType, 
		   int ruleId, 
		   int cellId)
{
  cell[x][y][z].set (aCellType, ruleId, cellId);
  numberOfCells++;
  if (x>maxX) maxX = x;
  if (y>maxY) maxY = y;
  if (z>maxZ) maxZ = z;
  if (x<minX) minX = x;
  if (y<minY) minY = y;
  if (z<minZ) minZ = z;
}

typeOfCell Pattern::getTypeOfCell (int x, int y, int z)
{
  return cell[x][y][z].value ();
}

int Pattern::getRuleId (int x, int y, int z)
{
  return cell[x][y][z].getAppliedRule ();
}

int Pattern::cellsNumber ()
{
  return numberOfCells;
}

void Pattern::debug ()
{
  int i, j, k;

  if (numberOfCells == 0) {
    printf ("empty pattern\n");
    return;
  }
  
  for (k=minZ; k<=maxZ; k++) {
    for (j=minY; j<=maxY; j++) {
      for (i=minX; i<=maxX; i++) {
	if (cell[i][j][k].value () == RED) {
	  if (cell[i][j][k].getAppliedRule () == initRuleId) {
	    printf ("R*");
	  }
	  else {
	    printf ("R ");
	  }
	}
	else if (cell[i][j][k].value () == YELLOW) {
	  if (cell[i][j][k].getAppliedRule () == initRuleId) {
	    printf ("Y*");
	  }
	  else {
	    printf ("Y ");
	  }
	}
	else {
	  printf ("- ");
	}
      }
      printf ("\n");
    }
    printf ("\n");
  }
  printf ("\n");
  
}

int Pattern::getMinX ()
{
  return minX;
}

int Pattern::getMaxX ()
{
  return maxX;
}

int Pattern::getMinY ()
{
  return minY;
}

int Pattern::getMaxY ()
{
  return maxY;
}

int Pattern::getMinZ ()
{
  return minZ;
}

int Pattern::getMaxZ ()
{
  return maxZ;
}

void Pattern::setEncoding ()
{
  int i;
  int lz;

  isEncoded = 1;
  lz = maxZ-minZ+1;

  /* on alloue les plans, le codage se fait dans le constructeur */

  plans = new PatternPlanEncoding*[lz];

  for (i=0; i<lz; i++) {
    plans[i] = new PatternPlanEncoding (this, i+minZ, simuType);
  }

}

void Pattern::freeEncoding ()
{
  int i;
  int lz;

  isEncoded = 0;
  lz = maxZ-minZ+1;

  /* et on libere les plans */

  for (i=0; i<lz; i++) {
    delete plans[i];
  }
  delete[] plans;

}

PatternPlanEncoding* Pattern::getPlan (int z)
{
  if ((z>maxZ-minZ) || (z<0)) {
    return NULL;
  }
  else {
    return plans[z];
  }
}

float Pattern::matchWithOtherPattern (Pattern* aPattern, 
				      correlationType* correlation, 
				      int* angle, 
				      int* dx,
				      int* dy,
				      int* dz)
{
  int d1, d2, d3;
  int i, j, k, ang;
  int mini;
  int bestd1, bestd2, bestd3;
  float correlationFitness, D1, D2, D3, N1, N2;
  float cf1, cf2;
  
  mini = cellsNumber()+aPattern->cellsNumber ();
  bestd1 = 0; /* needed by g++ -O2 */
  bestd2 = 0; /* needed by g++ -O2 */
  bestd3 = 0; /* needed by g++ -O2 */

#define matchingDepth 2

  for (ang=0; ang<4; ang++) {
    for (i=-matchingDepth; i<=matchingDepth; i++) {
      for (j=-matchingDepth; j<=matchingDepth; j++) {
	for (k=-matchingDepth; k<=matchingDepth; k++) {
	  matchPattern (aPattern, ang, i, j, k, &d1, &d2, &d3);
	  /* printf ("Resultat: %d d1=%d d2=%d d3=%d (a=%d dx=%d dy=%d dz=%d)\n",
	     d1+d2+d3, d1, d2, d3, ang, i, j, k); */
	  if (d1+d2+d3<mini) {
	    mini = d1+d2+d3;
	    bestd1=d1;
	    bestd2=d2;
	    bestd3=d3;
	    *angle=ang;
	    *dx = i;
	    *dy = j;
	    *dz = k;
	  }
	}
      }
    }    
  }

  /* calcul du type det de l'indice de correlation */

  D1 = (float)bestd1;
  D2 = (float)bestd2;
  D3 = (float)bestd3;
  N1 = (float)cellsNumber();
  N2 = (float)aPattern->cellsNumber ();

  if (bestd1+bestd2+bestd3 == 0) {
    *correlation = matches;
    correlationFitness = 1.0;
  }
  else if (bestd1+bestd3 == 0) {
    *correlation = includes;
    correlationFitness = 1.0;
  }
  else if (bestd2+bestd3 == 0) {
    *correlation = isIncluded;
    correlationFitness = 1.0;
  }
  else {
    *correlation = different;
    if ((D1 !=0) || (D2 !=0)) {
      cf1 = 1-(D1)/(D2+D1);
      cf2 = 1-(D2)/(D1+D2);
      if (cf1 > cf2) {
	correlationFitness = cf1;
      }
      else {
	correlationFitness = cf2;
      }
    }
    else {
      correlationFitness = 1-D3/(N1+N2);
    }
    correlationFitness = 1;
    correlationFitness = 1-(D1+D2+D3)/(N1+N2);
    
    correlationFitness = (float)(exp(log((double)correlationFitness)*100/(double)(N1+N2)));
  }
  
  /* des essais infructueux
     
     correlationFitness = 1-(D1+D2+D3)/(N1+N2);
     correlationFitness = 1-(D2*(D1+D3)+D1*(D2+D3))/(N1*N2);
     correlationFitness = correlationFitness*correlationFitness;
     correlationFitness = 1-D3/(N1+N2-D1-D2);
     correlationFitness = 1-(D1+D3)*(D2+D3)/(N1*N2);
     correlationFitness = 1-(D1+D2+D3)*D1*D2/(N1+N2)/(N1*N2);
     correlationFitness = (1-(D1+D2)/(N1+N2))*(1-D3/(N1+N2)); */
  
  /* printf (" mini=%d d1=%d d2=%d d3=%d - (n1=%d n2=%d) - correlation: %f \n",
     mini, bestd1, bestd2, bestd3, cellsNumber (), aPattern->cellsNumber (),
     correlationFitness); */
  
  return correlationFitness;
}

void Pattern::matchPattern (Pattern* aPattern,
			    int angle, int dx, int dy, int dz,
			    int* d1, int* d2, int* d3)
{
  int i;
  int max;
  int l1, l2;
  int ad1, ad2, ad3;
  PatternPlanEncoding* p1;
  PatternPlanEncoding* p2;
  
  l1 = maxZ-minZ+1;
  l2 = aPattern->getMaxZ()-aPattern->getMinZ()+1;
  ad1 = 0;
  ad2 = 0;
  ad3 = 0;
  *d1 = 0;
  *d2 = 0;
  *d3 = 0;

  if (dz<0) {
    
    if (l1>l2-dz) {
      max = l1;
    }
    else {
      max = l2-dz;
    }

    for (i=0; i<max; i++) {
      p1 = getPlan (i);
      p2 = aPattern->getPlan (i+dz);
      if ((p1 == NULL) && (p2 == NULL)) {
	ad1 = 0;
	ad2 = 0;
	ad3 = 0;
      }
      else if (p1 == NULL) {
	ad1 = p2->cellsNumber ();
	ad2 = 0;
	ad3 = 0;
      }
      else if (p2 == NULL) {
	ad1 = 0;
	ad2 = p1->cellsNumber ();
	ad3 = 0;
      }
      else {
	p1->matchWithOtherPlan (p2, angle, dx, dy, &ad1, &ad2, &ad3);
      }
      *d1 = *d1+ad1;
      *d2 = *d2+ad2;
      *d3 = *d3+ad3;
    }

  }
  
  else { /* dz>=0 */
    
    if (l2>l1+dz) {
      max = l2;
    }
    else {
      max = l1+dz;
    }

    for (i=0; i<max; i++) {
      p1 = getPlan (i-dz);
      p2 = aPattern->getPlan (i);
      if ((p1 == NULL) && (p2 == NULL)) {
	ad1 = 0;
	ad2 = 0;
	ad3 = 0;
      }
      else if (p1 == NULL) {
	ad1 = p2->cellsNumber ();
	ad2 = 0;
	ad3 = 0;
      }
      else if (p2 == NULL) {
	ad1 = 0;
	ad2 = p1->cellsNumber ();
	ad3 = 0;
      }
      else {
	p1->matchWithOtherPlan (p2, angle, dx, dy, &ad1, &ad2, &ad3);
      }
      *d1 = *d1+ad1;
      *d2 = *d2+ad2;
      *d3 = *d3+ad3;
    }

  }
    
}
			      
void Pattern::initWithPattern (Pattern* aPattern)
{
  int i, j, k;
  
  for (i=0; i<maxSize; i++) {
    for (j=0; j<maxSize; j++) {
      for (k=0; k<maxSize; k++) {
	cell[i][j][k].set (EMPTY, -1, -1);
      }
    }
  }
  
  minX = maxSize;
  minY = maxSize;
  minZ = maxSize;
  maxX = -1;
  maxY = -1;
  maxZ = -1;
  numberOfCells = 0;
  isEncoded = 0;
  classIndex = -1;

  for (i=aPattern->getMinX(); i<=aPattern->getMaxX (); i++) {
    for (j=aPattern->getMinY(); j<=aPattern->getMaxY (); j++) {
      for (k=aPattern->getMinZ(); k<=aPattern->getMaxZ (); k++) {
	set (i, j, k, aPattern->getTypeOfCell (i, j, k), 0, 0);
      }
    }
  }
  
}

void Pattern::mergeWithPattern (Pattern* aPattern,
				int angle, int dx, int dy, int dz)
{
  int i, j, k;
  typeOfCell cellType;

  printf ("Merging base pattern:\n");
  debug ();
  printf ("and pattern:\n");
  aPattern->debug ();
  printf ("angle: %d dx: %d dy: %d dz: %d\n", angle, dx, dy, dz);
  
  for (i=0; i<=aPattern->getMaxX ()-aPattern->getMinX(); i++) {
    for (j=0; j<=aPattern->getMaxY ()-aPattern->getMinY(); j++) {
      for (k=0; k<=aPattern->getMaxZ ()-aPattern->getMinZ(); k++) {
	
	if (angle == 0) {
	  cellType = aPattern->getTypeOfCell (i+aPattern->getMinX(),
					      j+aPattern->getMinY(),
					      k+aPattern->getMinZ());
	  if (cellType != EMPTY) {
	    set (minX+i-dx, minY+j-dy, minZ+k-dz, cellType, 0, 0);
	  }
	}
      }
    }
  }
  
  printf ("Resultat:\n");
  debug ();
 
  printf ("Pas fini d'etre implemente\n");
  exit (-1);
}

/*****************************/
/* Objet PatternPlanEncoding */
/*****************************/

PatternPlanEncoding::PatternPlanEncoding (Pattern* aPattern, 
					  int aZ, 
					  SimulationType aSimuType)
{
  int i, j;
  int lx, ly;
  int cellType;

  pattern = aPattern;
  z = aZ;
  simuType = aSimuType;
  numberOfCells = 0;

  if ((simuType == cubic) /* coding for cubic structures */
      || (simuType == hexa)) { /* coding for hexa structures */
    
    /* pour le moment, je n'ai pas le temps de proceder au codage des
       patterns en hexa en tenant compte des 6 rotations, alors, je code
       parel qu'en cubique */

    lx = maxX()-minX()+1;
    ly = maxY()-minY()+1;
    
    /* reservation de la place-memoire */

    cells = new int*[4];
    types = new int*[4];

    cells[0] = new int[ly];
    for (i=0; i<ly; i++) cells[0][i] = 0;
    cells[1] = new int[ly];
    for (i=0; i<ly; i++) cells[1][i] = 0;
    cells[2] = new int[lx];
    for (i=0; i<lx; i++) cells[2][i] = 0;
    cells[3] = new int[lx];
    for (i=0; i<lx; i++) cells[3][i] = 0;

    types[0] = new int[ly];
    for (i=0; i<ly; i++) types[0][i] = 0;
    types[1] = new int[ly];
    for (i=0; i<ly; i++) types[1][i] = 0;
    types[2] = new int[lx];
    for (i=0; i<lx; i++) types[2][i] = 0;
    types[3] = new int[lx];
    for (i=0; i<lx; i++) types[3][i] = 0;

    /* codage proprement dit */

    for (i=0; i<lx; i++) {
      for (j=0; j<ly; j++) {
	cellType = pattern->getTypeOfCell (i+minX(), j+minY(), z);
	/* NOTA: pour le moment, dans les patterns, RED=BLUE */
	if ((cellType == RED) || (cellType == BLUE)) {
	  cells[0][j]      += (1 << i);
	  cells[1][ly-j-1] += (1 << (lx-i-1));
	  cells[2][lx-i-1] += (1 << j);
	  cells[3][i]      += (1 << (ly-j-1));
	  types[0][j]      += (1 << i);
	  types[1][ly-j-1] += (1 << (lx-i-1));
	  types[2][lx-i-1] += (1 << j);
	  types[3][i]      += (1 << (ly-j-1));
	  numberOfCells++;
	}
	else if (cellType == YELLOW) {
	  cells[0][j]      += (1 << i);
	  cells[1][ly-j-1] += (1 << (lx-i-1));
	  cells[2][lx-i-1] += (1 << j);
	  cells[3][i]      += (1 << (ly-j-1));
	  numberOfCells++;
	}
      }
    }
    
    /* (pas facile a lire tout ca, hein ???
       methode: se concentrer un peu, et ca vient tout seul) */
    
  }
  else { /* coding for hexa structures */ 
    printf ("Sorry: PatternPlanEncoding not implemented in hexa\n");
  }
}

PatternPlanEncoding::~PatternPlanEncoding ()
{
  if ((simuType == cubic) /* free for cubic structures */
      || (simuType == hexa)) { /* free for hexa structures */
    delete[] cells[0];
    delete[] cells[1];
    delete[] cells[2];
    delete[] cells[3];
    delete[] types[0];
    delete[] types[1];
    delete[] types[2];
    delete[] types[3];
    delete[] cells;
    delete[] types;
  }
  else { /* free for hexa structures */
    printf ("Sorry: PatternPlanEncoding not implemented in hexa\n");
  }
}

SimulationType PatternPlanEncoding::getSimuType ()
{
  return simuType;
}

int PatternPlanEncoding::getCells (int angle, int index)
{
  if ((simuType == cubic) /* getCells for cubic structures */
      || (simuType == hexa)) /* pareil */
      { 

	  if (index<0)
	      {
		  return 0;
	      }
	  else
	      {
		  if ((angle == 0) || (angle == 1))
		      {
			  if (index > (maxY()-minY()))
			      return 0;
			  else
			      return cells[angle][index];
		      }
		  else /* angle = 2 ou angle = 3*/
		      { 
			  if (index > (maxX()-minX()))
			      return 0;
			  else
			      return cells[angle][index];
		      }
	      }
      }
  
  else
      { /* getCells for hexa structures */
      }
  
  return -1;
}

int PatternPlanEncoding::getTypes (int angle, int index)
{
  if ((simuType == cubic) /* getTypes for cubic structures */
      || (simuType == hexa)) /* pareil */

      {
	  if (index<0) 
	      return 0;
	  
	  else
	      {
		  if ((angle == 0) || (angle == 1))
		      {
			  if (index > (maxY()-minY()))
			      return 0;
			  else 
			      return types[angle][index];
		      }
		  else /* angle = 2 ou angle = 3*/
		      { 
			  if (index > (maxX()-minX()))
			      return 0;
			  else
			      return types[angle][index];
		      }
	      }
      }

  else /* getTypes for hexa structures */
      { 
      }
  
  return -1;
}

int PatternPlanEncoding::nbits (int anInt)
{
  int i;
  int returned;
  int current;

  returned = 0;
  current= anInt;
  for (i=0; i<maxLength*4; i++) {
    if (current & 1) returned++;
    current = current >> 1;
  }
  
  return returned; 
}

int PatternPlanEncoding::cellsNumber ()
{
  return numberOfCells;
}

void PatternPlanEncoding::matchWithOtherPlan (PatternPlanEncoding* aPlan, 
					      int angle,
					      int dx,
					      int dy,
					      int* d1,
					      int* d2,
					      int* d3)
/* les trois dernieres variables sont retournees et correspondent:
   d1: au nombre de cellules "en moins" dans le patternPlan courant
   d2: au nombre de cellules "en moins" dans le patternPlan passe en parametre
   d3: au nombre de cellules qui matchent bien, mais pour lesquelles on
   a un probleme de couleur */

/* P1: patternPlan courant
   P2: patternPlan passe en parametre
   c1 = cells (P1) t1 = types (P1)
   c2 = cells (P2) t2 = types (P2)
   d1 = nbits (~c1 & c2)
   d2 = nbits (c1 & ~c2)
   d3 = nbits (c1 & c2 & [(~t1 & t2) | (t1 & ~t2)]) */
{
  int c1, t1;
  int c2, t2;
  int max;
  int l1, l2;
  int j;

  if (simuType != aPlan->getSimuType ()) {
    printf ("ERROR: try to match two differents structures types (cubic and hexa)\n");
  }
  else {
    if ((simuType == cubic)  /* pattern matching for cubic structures */
	|| (simuType == hexa)) { /* pareil */

      /* un peu de debug pour comprendre */

      /* printf ("\nOn compare avec: angle=%d dx=%d dy=%d\n", angle, dx, dy);
	 printf ("le plan:\n");
	 debug (0);
	 printf ("et le plan:\n");
	 aPlan->debug (angle); */

      /* on y va maintenant... (c'est une partie dure, la) */
      
      *d1 = 0;
      *d2 = 0;
      *d3 = 0;
      
      l1 = maxY()-minY()+1;
      if ((angle == 0) || (angle == 1)) {
	l2 = aPlan->maxY()-aPlan->minY()+1;
      }
      else {
	l2 = aPlan->maxX()-aPlan->minX()+1;
      }
      
      /* NOTE: l'angle est applique au patternPlan passe en parametre,
	 les decalages sont appliques aux patternPlans en fonction de leur
	 signe */
      
      if (dy<0) {
	if (l1 > l2-dy) {
	  max = l1;
	}
	else {
	  max = l2-dy;
	}
	for (j=0; j<max ; j++) {
	  c1 = getCells (0, j);
	  t1 = getTypes (0, j);
	  c2 = aPlan->getCells (angle, j+dy);
	  t2 = aPlan->getTypes (angle, j+dy);
	  if (dx<0) {
	    c1 = c1 << (-dx);
	    t1 = t1 << (-dx);
	  }
	  else {
	    c2 = c2 << dx;
	    t2 = t2 << dx;
	  }

	  *d1 += nbits (~(c1) & c2);
	  *d2 += nbits (c1 & ~(c2));
	  *d3 += nbits (c1 & c2 & (((~t1) & t2) | (t1 & (~t2))));

	  /* printf ("j=%d (c1=%d t1=%d c2=%d t2=%d) / d1=%d d2=%d d3=%d\n",
	     j, c1, t1, c2, t2, 
	     nbits (~(c1) & c2), 
	     nbits (c1 & ~(c2)),
	     nbits (c1 & c2 & (((~t1) & t2) | (t1 & (~t2))))); */
	  
	}
      }
      
      else {
	if (l1+dy > l2) {
	  max = l1+dy;
	}
	else {
	  max = l2;
	}
	for (j=0; j<max ; j++) {
	  c1 = getCells (0, j-dy);
	  t1 = getTypes (0, j-dy);
	  c2 = aPlan->getCells (angle, j);
	  t2 = aPlan->getTypes (angle, j);
	  if (dx<0) {
	    c1 = c1 << (-dx);
	    t1 = t1 << (-dx);
	  }
	  else {
	    c2 = c2 << dx;
	    t2 = t2 << dx;
	  }
	  
	  *d1 += nbits (~(c1) & c2);
	  *d2 += nbits (c1 & ~(c2));
	  *d3 += nbits (c1 & c2 & (((~t1) & t2) | (t1 & (~t2))));

	  /* printf ("j=%d (c1=%d t1=%d c2=%d t2=%d) / d1=%d d2=%d d3=%d\n",
	     j, c1, t1, c2, t2, 
	     nbits (~(c1) & c2), 
	     nbits (c1 & ~(c2)),
	     nbits (c1 & c2 & (((~t1) & t2) | (t1 & (~t2))))); */
	  
	}
      }
      
      /* printf ("On trouve: d1=%d, d2=%d, d3=%d\n", *d1, *d2, *d3); */

    }
    else {                   /* pattern matching for hexa structures */
      printf ("Sorry: PatternPlanEncoding not implemented in hexa\n");
    }
  }
}

void PatternPlanEncoding::debug (int angle)
{
  int i, j;
  int lx, ly;

  if ((simuType == cubic)  /* debugging for cubic structures */
      || (simuType == hexa)) /* pareil */
      { 
    
	  lx = maxX()-minX()+1;
	  ly = maxY()-minY()+1;
	  switch (angle)
	      {
	      case 0:
		  printf ("Debug PatternPlanEncoding for rotation: 0\n");
		  for (j=0; j<ly; j++) 
		      printf ("cells %d types %d\n", cells[0][j], types[0][j]);
		  break;
		  
	      case 1:
		  printf ("Debug PatternPlanEncoding for rotation: 1\n");
		  for (j=0; j<ly; j++)
		      printf ("cells %d types %d\n", cells[1][j], types[1][j]);
		  break;
		  
	      case 2:
		  printf ("Debug PatternPlanEncoding for rotation: 2\n");
		  for (i=0; i<lx; i++) 
		      printf ("cells %d types %d\n", cells[2][i], types[2][i]);
		  break;
		  
	      default:
		  printf ("Debug PatternPlanEncoding for rotation: 3\n");
		  for (i=0; i<lx; i++) 
		      printf ("cells %d types %d\n", cells[3][i], types[3][i]);
		  
	      }
      }
  
  else
      { /* debugging for hexa structures */
	  printf ("Sorry: PatternPlanEncoding not implemented in hexa\n");
      }

}

int PatternPlanEncoding::minX ()
{
  return pattern->getMinX ();
}

int PatternPlanEncoding::maxX ()
{
  return pattern->getMaxX ();
}
     
int PatternPlanEncoding::minY ()
{
  return pattern->getMinY ();
}
     
int PatternPlanEncoding::maxY ()
{
  return pattern->getMaxY ();
}
     
int PatternPlanEncoding::minZ ()
{
  return pattern->getMinZ ();
}
     
int PatternPlanEncoding::maxZ ()
{
  return pattern->getMaxZ ();
}