Logo Search packages:      
Sourcecode: randim version File versions  Download package

drawprep.c

/* drawprep.c: code for specifying transformations, and computing
   scaling and relative probabilities for randim. */

/* 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.  Please see the file "COPYING"
   for details.  If you have not received it along with this program,
   please write to the Free Software Foundation, Inc., 675 Mass. Ave.,
   Cambridge, MA 02139, USA. */

#include "randim.h"

/**********************************************************************/
/* transformation specification */

/*This function generates a set of affine transformations based on the
  current parameters.*/
void
define(int trans_type)
{
  double    angle,angle2,r,r2;      /*for rotational type*/
  int       i,j,k;                  /*indices*/

  sqrand(seed);         /*Re-seed and increment each time so*/
  srandom(seed);        /*runs can be recovered.*/
  seed += 2;

  switch (trans_type) {
    
  case 1:   /*Fully random.*/
    for (i=0;i<n_trans;i++)
      for (j=0;j < nparams(dimension);j++)
      a[i][j] = dr() - 0.5;
    break;

  case 2: /*Fully symmetric.*/
    /*(Alternate = same entries all signs reversed.)*/
    for (i=0;i<n_trans;i+=2)
      for (j=0;j < nparams(dimension);j++) {
      a[i][j] = dr() - 0.5;
      a[i+1][j] = -a[i][j];
        }
    break;
    
  case 3:   /*Randomized symmetric.*/
    /*(Groups of 4, last 3 random +/- entries of first.)*/
    for (i=0;i<n_trans;i+=4)
      for (j=0;j < nparams(dimension);j++) {
      a[i][j] = dr() - 0.5;
      for (k=1;k<4;k++) {
        a[i+k][j] = (dr() < 0.5) ?
          a[i][j] : -a[i][j];
          }
        }
    break;
    
  case 4: /*Rotational, constant radius.*/
    /*(Rotation through a random angle [0,2pi), no scaling.)*/
    if (dimension == 2) {
      for (i=0;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr();
      a[i][0] = a[i][4] = cos(angle);
      a[i][3] = sin(angle);
      a[i][1] = -a[i][3];
      a[i][2] = dr() - 0.5;
      a[i][5] = dr() - 0.5;
        }
      }
    else {  /*dimension = 3*/
      /* (Here choose two angles (azimuth, elevation).)*/
      for (i=0;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr(); angle2 = 2.0 *
                             M_PI * dr();
      a[i][0] = cos(angle);
      a[i][1] = -sin(angle) * cos(angle2);
      a[i][2] = -sin(angle) * sin(angle2);
      a[i][4] = sin(angle) * cos(angle2);
      a[i][5] = cos(angle) * cos(angle2);
      a[i][6] = cos(angle) * sin(angle2);
      a[i][8] = sin(angle) * sin(angle2);
      a[i][9] = cos(angle) * sin(angle2);
      a[i][10]= cos(angle2);
      a[i][3] = dr(); a[i][7] = dr() - 0.5;
      a[i][11] = dr() - 0.5;
        }
      }
    break;
    
  case 5:   /*Rotational, randomized radius.*/
    /*(Same as above but scale by factor in [0,1).)*/
    if (dimension == 2) {
      for (i=0;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr();
      r = dr();
      a[i][0] = a[i][4] = r * cos(angle);
      a[i][3] = r * sin(angle);
      a[i][1] = -a[i][3];
      a[i][2] = dr() - 0.5;
      a[i][5] = dr() - 0.5;
        }
      }
    else {  /*dimension = 3*/
      for (i=0;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr(); angle2 = 2.0 *
                             M_PI * dr();
      r = dr(); r2 = dr();
      a[i][0] = r * cos(angle);
      a[i][1] = r * r2 * -sin(angle) * cos(angle2);
      a[i][2] = r * r2 * -sin(angle) * sin(angle2);
      a[i][4] = r * r2 * sin(angle) * cos(angle2);
      a[i][5] = r * r2 * cos(angle) * cos(angle2);
      a[i][6] = r * r2 * cos(angle) * sin(angle2);
      a[i][8] = r * r2 * sin(angle) * sin(angle2);
      a[i][9] = r * r2 * cos(angle) * sin(angle2);
      a[i][10]= r2 * cos(angle2);
      a[i][3] = dr(); a[i][7] = dr() - 0.5;
      a[i][11] = dr() - 0.5;
        }
      }
    break;
    
  case 6:   /*Hybrid.*/
    for (i=0;i<n_trans - hfactor;i++)
      for (j=0;j<nparams(dimension);j++)
      a[i][j] = dr() - 0.5;
    if (dimension == 2) {
      for (;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr();
      r = dr();
      a[i][0] = a[i][4] = r * cos(angle);
      a[i][3] = r * sin(angle);
      a[i][1] = -a[i][3];
      a[i][2] = dr() - 0.5;
      a[i][5] = dr() - 0.5;
        }

      }
    else {  /*dimension = 3*/
      for (;i<n_trans;i++) {
      angle = 2.0 * M_PI * dr();
      angle2 = 2.0 * M_PI * dr();
      r = dr(); r2 = dr();
      a[i][0] = r * cos(angle);
      a[i][1] = r * r2 * -sin(angle) * cos(angle2);
      a[i][2] = r * r2 * -sin(angle) * sin(angle2);
      a[i][4] = r * r2 * sin(angle) * cos(angle2);
      a[i][5] = r * r2 * cos(angle) * cos(angle2);
      a[i][6] = r * r2 * cos(angle) * sin(angle2);
      a[i][8] = r * r2 * sin(angle) * sin(angle2);
      a[i][9] = r * r2 * cos(angle) * sin(angle2);
      a[i][10]= r2 * cos(angle2);
      a[i][3] = dr(); a[i][7] = dr();
      a[i][11] = dr();
        }
      }                       
    break;    
    
  case 7:  /*Fern - gotten from initialization.*/
    dimension = 2;      /*restricted*/
    for (i=0;i<4;i++)
      for (j=0;j<6;j++)
      a[i][j] = fern[i][j];
    break;
    
  default:printf("\n\nBad trans_type.");
    exit(1);
    break;

    }

    for (i=0;i<n_trans;i++) {
      if (dimension == 2) {
      for (j=0;j<6;j++)
        if (isnan(a[i][j])) printf("** isnan: %d, %d\n",i,j);
      }
      else {
      for (j=0;j<12;j++)
        if (isnan(a[i][j])) printf("** isnan: %d, %d\n",i,j);
      }
      }

}


/**********************************************************************/
/* drawing preparation */

/* scale(), given transformations, calculates scaling and offset.  It works by
  iterating through a large number of points on the attractor and keeping
  track of their coordinate maxima and minima.  These values are then used
  to compute the scaling.*/
void
scale(int n_pts)
{
  int             i;
  float     xmin,ymin,xmax,ymax,zmin,zmax;
  float           temp;


  /*Set probability of selection of transformation in accordance with
  magnitude of its determinant.*/
  calc_prob_array();

  /*Autoscale based on first n_pts points.
    Uses pointer incrementing for speed. Hence have arrangements
                        [1 2 3] [4]
  [1 2] [3]             [5 6 7] [8]
  [4 5] [6] for 2d,           [9 a b] [c] for 3d.     
*/
  
  gx = gy = gz = 0.1;
  xmax = ymax = xmin = ymin = 0.0;
  xscale = yscale = zscale = oxscale = oyscale = ozscale = 1.0;
  for (i=0;i<n_pts;i++) {
    apply_trans();
    xmax = max(xmax,vx);
    ymax = max(ymax,vy);
    zmax = max(zmax,gz);
    xmin = min(xmin,vx);
    ymin = min(ymin,vy);
    zmin = min(zmin,gz);
      }

  /* total hack */
  if (isnan(xmax) || isnan(xmin) || isnan(ymax) || isnan(ymin)
      || isnan(zmax) || isnan(zmin)) {
    xmin = -100, xmax = 100, ymin = -100, ymax = 100, zmin = -100, zmax = 100;
    printf("HACKED!\n");
  }

  /* note, z-coordinate scaled differently from x and y, to 1.0 */
  temp = (hsize / (xmax - xmin));
  oxscale = xscale = temp;
  temp = (vsize / (ymax - ymin));
  oyscale = yscale = temp;
  temp = (1.0 / (zmax - zmin));
  ozscale = zscale = temp;

  temp = -hsize2 - (xmin * xscale);
  oxoffset = xoffset = temp;        /*subtract so 0,0 will be center*/
  temp = -vsize2 - (ymin * yscale); /*of scaled image*/
  oyoffset = yoffset = temp;
  temp = -(zmin * zscale);
  ozoffset = zoffset = temp;
  zxmag = zymag = 1.0;              /*reset zoom parameters*/
  zoomf = 0;

  /*scalea little extra just in case*/
  oxscale /= 1.1, xscale /= 1.1, oyscale /= 1.1, yscale /= 1.1,
    ozscale /= 1.1, zscale /= 1.1;
}

/* calc_prob_array() sets up probability array (prob[]) according to the
   relative sizes of the determinants of the linear parts of the transform-
   ations.  Basically, the larger the determinant corresponding to a partic-
   ular index is, the more entries in the array will be that index.  Thus,
   when a random element is picked from the array during drawing, the probab-
   ilities are observed, and we avoid having to do a bunch of floating point
   if-tests every time. */
void
calc_prob_array()
{
  int       i,j,              /*indices*/
            index[TRANS_MAX+1];     /*prob indices where entry changes*/
  char            zeros=0, maxidx=0, minidx=0;  /*obvious purposes*/
  int       minval;
  float           det[TRANS_MAX+1], /*hold determinants*/
            total = 0.0;            /*for normalization*/

  for (i=0;i<n_trans;i++) {   /*calculate determinants*/
    if (dimension == 2)
      det[i+1] = a[i][0] * a[i][4] - a[i][1] * a[i][3];
    else    /*dimension = 3*/
      det[i+1] = a[i][0] * (a[i][5] * a[i][10] - a[i][6] * a[i][9]) -
             a[i][1] * (a[i][4] * a[i][10] - a[i][6] * a[i][8]) +
             a[i][2] * (a[i][4] * a[i][9] - a[i][5] * a[i][8]);
    det[i+1] = (det[i+1] > 0) ? det[i+1] : -det[i+1]; /*abs*/
    total += det[i+1];
    }

  for (i=1;i<n_trans+1;i++) { /*find #/0's and maximum det*/
    if (det[i] == 0)
      zeros++;
    else {
      maxidx = (det[i] < det[maxidx]) ? maxidx : i;
      minidx = (det[i] > det[minidx]) ? minidx : i;
      }
    }
  minval = det[minidx];


  index[0] = 0;
  for (i=1;i<n_trans;i++) {   /*find indices*/
    if (det[i] == 0)
      index[i] = index[i-1] + minval;     /*minimal*/
    else {
      index[i] = index[i-1] + (det[i] / total) * 32768.0;
      if (i == maxidx)
      index[i] -= zeros * minval;/*(arbitrarily take room for 0s from here)*/
      }
    }
  index[n_trans] = 32768;           /*so we fill to end of prob[] array*/

  for (i=1;i<n_trans+1;i++)   /*fill in prob[] according to index[]*/
    for (j=index[i-1];j<index[i];j++)
      prob[j] = i - 1;

  for (i=0;i<32768;i++) {
    if ((prob[i] > n_trans) || (prob[i] < 0))
      printf("calc_prob_array(): error %d\n",i);
    }

}

Generated by  Doxygen 1.6.0   Back to index