home *** CD-ROM | disk | FTP | other *** search
/ Encyclopedia of Graphics File Formats Companion / GFF_CD.ISO / formats / ttddd / spec / t3d_doc / igensurf.zoo / src / noise3.c < prev    next >
C/C++ Source or Header  |  1991-09-28  |  5KB  |  246 lines
  1. /* :ts=8 */
  2. /* Copyright (c) 1988 Regents of the University of California */
  3.  
  4. #ifndef lint
  5. static char SCCSid[] = "@(#)noise3.c 1.3 10/30/90 LBL";
  6. #endif
  7.  
  8. /*
  9.  *  noise3.c - noise functions for random textures.
  10.  *
  11.  *     Credit for the smooth algorithm goes to Ken Perlin.
  12.  *     (ref. SIGGRAPH Vol 19, No 3, pp 287-96)
  13.  *
  14.  *     4/15/86
  15.  *     5/19/88    Added fractal noise function
  16.  */
  17.  
  18.  
  19. #define  A        0
  20. #define  B        1
  21. #define  C        2
  22. #define  D        3
  23.  
  24. #define  rand3a(x,y,z)    frand(67*(x)+59*(y)+71*(z))
  25. #define  rand3b(x,y,z)    frand(73*(x)+79*(y)+83*(z))
  26. #define  rand3c(x,y,z)    frand(89*(x)+97*(y)+101*(z))
  27. #define  rand3d(x,y,z)    frand(103*(x)+107*(y)+109*(z))
  28.  
  29. #define  hermite(p0,p1,r0,r1,t)  (    p0*((2.0*t-3.0)*t*t+1.0) + \
  30.                     p1*(-2.0*t+3.0)*t*t + \
  31.                     r0*((t-2.0)*t+1.0)*t + \
  32.                     r1*(t-1.0)*t*t )
  33.     
  34. double  *noise3(), noise3coef(), Get_Argument(), frand();
  35.  
  36. static long  xlim[3][2];
  37. static double  xarg[3];
  38.  
  39. #define  EPSILON    .0001        /* error allowed in fractal */
  40.  
  41. #define  frand3(x,y,z)    frand(17*(x)+23*(y)+29*(z))
  42.  
  43. double  fnoise3();
  44.  
  45.  
  46. double
  47. l_noise3()            /* compute 3-dimensional noise function */
  48. {
  49.     return(noise3coef(D));
  50. }
  51.  
  52.  
  53. double
  54. l_noise3a()            /* compute x slope of noise function */
  55. {
  56.     return(noise3coef(A));
  57. }
  58.  
  59.  
  60. double
  61. l_noise3b()            /* compute y slope of noise function */
  62. {
  63.     return(noise3coef(B));
  64. }
  65.  
  66.  
  67. double
  68. l_noise3c()            /* compute z slope of noise function */
  69. {
  70.     return(noise3coef(C));
  71. }
  72.  
  73.  
  74. double
  75. l_fnoise3()            /* compute fractal noise function */
  76. {
  77.     double  x[3];
  78.  
  79.     x[0] = Get_Argument(1);
  80.     x[1] = Get_Argument(2);
  81.     x[2] = Get_Argument(3);
  82.  
  83.     return(fnoise3(x));
  84. }
  85.  
  86.  
  87. static double
  88. noise3coef(coef)        /* return coefficient of noise function */
  89. int  coef;
  90. {
  91.     double  x[3];
  92.  
  93.     x[0] = Get_Argument(1);
  94.     x[1] = Get_Argument(2);
  95.     x[2] = Get_Argument(3);
  96.  
  97.     return(noise3(x)[coef]);
  98. }
  99.  
  100.  
  101. double *
  102. noise3(xnew)            /* compute the noise function */
  103. register double  xnew[3];
  104. {
  105.     extern double  floor();
  106.     static double  x[3] = {-100000.0, -100000.0, -100000.0};
  107.     static double  f[4];
  108.  
  109.     if (x[0]==xnew[0] && x[1]==xnew[1] && x[2]==xnew[2])
  110.         return(f);
  111.     x[0] = xnew[0]; x[1] = xnew[1]; x[2] = xnew[2];
  112.     xlim[0][0] = floor(x[0]); xlim[0][1] = xlim[0][0] + 1;
  113.     xlim[1][0] = floor(x[1]); xlim[1][1] = xlim[1][0] + 1;
  114.     xlim[2][0] = floor(x[2]); xlim[2][1] = xlim[2][0] + 1;
  115.     xarg[0] = x[0] - xlim[0][0];
  116.     xarg[1] = x[1] - xlim[1][0];
  117.     xarg[2] = x[2] - xlim[2][0];
  118.     interpolate(f, 0, 3);
  119.     return(f);
  120. }
  121.  
  122.  
  123. static
  124. interpolate(f, i, n)
  125. double  f[4];
  126. register int  i, n;
  127. {
  128.     double  f0[4], f1[4];
  129.  
  130.     if (n == 0) {
  131.         f[A] = rand3a(xlim[0][i&1],xlim[1][i>>1&1],xlim[2][i>>2]);
  132.         f[B] = rand3b(xlim[0][i&1],xlim[1][i>>1&1],xlim[2][i>>2]);
  133.         f[C] = rand3c(xlim[0][i&1],xlim[1][i>>1&1],xlim[2][i>>2]);
  134.         f[D] = rand3d(xlim[0][i&1],xlim[1][i>>1&1],xlim[2][i>>2]);
  135.     } else {
  136.         n--;
  137.         interpolate(f0, i, n);
  138.         interpolate(f1, i | 1<<n, n);
  139.         f[A] = (1.0-xarg[n])*f0[A] + xarg[n]*f1[A];
  140.         f[B] = (1.0-xarg[n])*f0[B] + xarg[n]*f1[B];
  141.         f[C] = (1.0-xarg[n])*f0[C] + xarg[n]*f1[C];
  142.         f[D] = hermite(f0[D], f1[D], f0[n], f1[n], xarg[n]);
  143.     }
  144. }
  145.  
  146.  
  147. double
  148. frand(s)            /* get random number from seed */
  149. register long  s;
  150. {
  151.     s = s<<13 ^ s;
  152.     return(1.0-((s*(s*s*15731+789221)+1376312589)&0x7fffffff)/1073741824.0);
  153. }
  154.  
  155.  
  156. double
  157. l_hermite()            /* library call for hermite interpolation */
  158. {
  159.     double  t;
  160.     
  161.     t = Get_Argument(5);
  162.     return( hermite(Get_Argument(1), Get_Argument(2),
  163.             Get_Argument(3), Get_Argument(4), t) );
  164. }
  165.  
  166.  
  167. double
  168. fnoise3(p)            /* compute fractal noise function */
  169. double  p[3];
  170. {
  171.     double  floor();
  172.     long  t[3], v[3], beg[3], s;
  173.     double  fval[8], fc;
  174.     int  branch;
  175.     register int  i, j;
  176.                         /* get starting cube */
  177.     s = (long)(1.0/EPSILON);
  178.     for (i = 0; i < 3; i++) {
  179.         t[i] = s*p[i];
  180.         beg[i] = s*floor(p[i]);
  181.     }
  182.     for (j = 0; j < 8; j++) {
  183.         for (i = 0; i < 3; i++) {
  184.             v[i] = beg[i];
  185.             if (j & 1<<i)
  186.                 v[i] += s;
  187.         }
  188.         fval[j] = frand3(v[0],v[1],v[2]);
  189.     }
  190.                         /* compute fractal */
  191.     for ( ; ; ) {
  192.         s >>= 1;
  193.         branch = 0;
  194.         for (i = 0; i < 3; i++) {    /* do center */
  195.             v[i] = beg[i] + s;
  196.             if (t[i] > v[i]) {
  197.                 branch |= 1<<i;
  198.             }
  199.         }
  200.         fc = 0.0;
  201.         for (j = 0; j < 8; j++)
  202.             fc += fval[j];
  203.         fc *= 0.125;
  204.         if (s < 1)
  205.             return(fc);        /* close enough */
  206.         fc += s*EPSILON*frand3(v[0],v[1],v[2]);
  207.         fval[~branch & 7] = fc;
  208.         for (i = 0; i < 3; i++) {    /* do faces */
  209.             if (branch & 1<<i)
  210.                 v[i] += s;
  211.             else
  212.                 v[i] -= s;
  213.             fc = 0.0;
  214.             for (j = 0; j < 8; j++)
  215.                 if (~(j^branch) & 1<<i)
  216.                     fc += fval[j];
  217.             fc = 0.25*fc + s*EPSILON*frand3(v[0],v[1],v[2]);
  218.             fval[~(branch^1<<i) & 7] = fc;
  219.             v[i] = beg[i] + s;
  220.         }
  221.         for (i = 0; i < 3; i++) {    /* do edges */
  222.             j = (i+1)%3;
  223.             if (branch & 1<<j)
  224.                 v[j] += s;
  225.             else
  226.                 v[j] -= s;
  227.             j = (i+2)%3;
  228.             if (branch & 1<<j)
  229.                 v[j] += s;
  230.             else
  231.                 v[j] -= s;
  232.             fc = fval[branch & ~(1<<i)];
  233.             fc += fval[branch | 1<<i];
  234.             fc = 0.5*fc + s*EPSILON*frand3(v[0],v[1],v[2]);
  235.             fval[branch^1<<i] = fc;
  236.             j = (i+1)%3;
  237.             v[j] = beg[j] + s;
  238.             j = (i+2)%3;
  239.             v[j] = beg[j] + s;
  240.         }
  241.         for (i = 0; i < 3; i++)        /* new cube */
  242.             if (branch & 1<<i)
  243.                 beg[i] += s;
  244.     }
  245. }
  246.