ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/noise3.c
(Generate patch)

Comparing ray/src/rt/noise3.c (file contents):
Revision 1.3 by greg, Tue Oct 30 21:03:49 1990 UTC vs.
Revision 2.12 by greg, Tue Feb 22 16:45:12 2011 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1988 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   *  noise3.c - noise functions for random textures.
6   *
7   *     Credit for the smooth algorithm goes to Ken Perlin.
8   *     (ref. SIGGRAPH Vol 19, No 3, pp 287-96)
12 *
13 *     4/15/86
14 *     5/19/88  Added fractal noise function
9   */
10  
11 + #include "copyright.h"
12  
13 + #include  "ray.h"
14 + #include  "func.h"
15 +
16   #define  A              0
17   #define  B              1
18   #define  C              2
# Line 25 | Line 23 | static char SCCSid[] = "$SunId$ LBL";
23   #define  rand3c(x,y,z)  frand(89*(x)+97*(y)+101*(z))
24   #define  rand3d(x,y,z)  frand(103*(x)+107*(y)+109*(z))
25  
26 < #define  hermite(p0,p1,r0,r1,t)  (      p0*((2.0*t-3.0)*t*t+1.0) + \
27 <                                        p1*(-2.0*t+3.0)*t*t + \
28 <                                        r0*((t-2.0)*t+1.0)*t + \
29 <                                        r1*(t-1.0)*t*t )
26 > #define  hpoly1(t)      ((2.0*t-3.0)*t*t+1.0)
27 > #define  hpoly2(t)      (-2.0*t+3.0)*t*t
28 > #define  hpoly3(t)      ((t-2.0)*t+1.0)*t
29 > #define  hpoly4(t)      (t-1.0)*t*t
30  
31 < double  *noise3(), noise3coef(), argument(), frand();
31 > #define  hermite(p0,p1,r0,r1,t)  (      p0*hpoly1(t) + \
32 >                                        p1*hpoly2(t) + \
33 >                                        r0*hpoly3(t) + \
34 >                                        r1*hpoly4(t) )
35  
36 + static char  noise_name[4][8] = {"noise3x", "noise3y", "noise3z", "noise3"};
37 + static char  fnoise_name[] = "fnoise3";
38 + static char  hermite_name[] = "hermite";
39 +
40   static long  xlim[3][2];
41   static double  xarg[3];
42  
43 < #define  EPSILON        .0001           /* error allowed in fractal */
43 > #define  EPSILON        .001            /* error allowed in fractal */
44  
45   #define  frand3(x,y,z)  frand(17*(x)+23*(y)+29*(z))
46  
47 < double  fnoise3();
47 > static double l_noise3(char  *nam);
48 > static double l_hermite(char *nm);
49 > static double * noise3(double  xnew[3]);
50 > static void interpolate(double  f[4], int  i, int  n);
51 > static double frand(long  s);
52 > static double fnoise3(double  p[3]);
53  
54  
55 < double
56 < l_noise3()                      /* compute 3-dimensional noise function */
55 > static double
56 > l_noise3(                       /* compute a noise function */
57 >        register char  *nam
58 > )
59   {
60 <        return(noise3coef(D));
49 < }
50 <
51 <
52 < double
53 < l_noise3a()                     /* compute x slope of noise function */
54 < {
55 <        return(noise3coef(A));
56 < }
57 <
58 <
59 < double
60 < l_noise3b()                     /* compute y slope of noise function */
61 < {
62 <        return(noise3coef(B));
63 < }
64 <
65 <
66 < double
67 < l_noise3c()                     /* compute z slope of noise function */
68 < {
69 <        return(noise3coef(C));
70 < }
71 <
72 <
73 < double
74 < l_fnoise3()                     /* compute fractal noise function */
75 < {
60 >        register int  i;
61          double  x[3];
62 <
62 >                                        /* get point */
63          x[0] = argument(1);
64          x[1] = argument(2);
65          x[2] = argument(3);
66 <
67 <        return(fnoise3(x));
66 >                                        /* make appropriate call */
67 >        if (nam == fnoise_name)
68 >                return(fnoise3(x));
69 >        i = 4;
70 >        while (i--)
71 >                if (nam == noise_name[i])
72 >                        return(noise3(x)[i]);
73 >        eputs(nam);
74 >        eputs(": called l_noise3!\n");
75 >        quit(1);
76 >        return 1; /* pro forma return */
77   }
78  
79  
80   static double
81 < noise3coef(coef)                /* return coefficient of noise function */
88 < int  coef;
81 > l_hermite(char *nm)             /* library call for hermite interpolation */
82   {
83 <        double  x[3];
83 >        double  t;
84 >        
85 >        t = argument(5);
86 >        return( hermite(argument(1), argument(2),
87 >                        argument(3), argument(4), t) );
88 > }
89  
92        x[0] = argument(1);
93        x[1] = argument(2);
94        x[2] = argument(3);
90  
91 <        return(noise3(x)[coef]);
91 > extern void
92 > setnoisefuncs(void)                     /* add noise functions to library */
93 > {
94 >        register int  i;
95 >
96 >        funset(hermite_name, 5, ':', l_hermite);
97 >        funset(fnoise_name, 3, ':', l_noise3);
98 >        i = 4;
99 >        while (i--)
100 >                funset(noise_name[i], 3, ':', l_noise3);
101   }
102  
103  
104 < double *
105 < noise3(xnew)                    /* compute the noise function */
106 < register double  xnew[3];
104 > static double *
105 > noise3(                 /* compute the noise function */
106 >        register double  xnew[3]
107 > )
108   {
104        extern double  floor();
109          static double  x[3] = {-100000.0, -100000.0, -100000.0};
110          static double  f[4];
111  
# Line 119 | Line 123 | register double  xnew[3];
123   }
124  
125  
126 < static
127 < interpolate(f, i, n)
128 < double  f[4];
129 < register int  i, n;
126 > static void
127 > interpolate(
128 >        double  f[4],
129 >        register int  i,
130 >        register int  n
131 > )
132   {
133 <        double  f0[4], f1[4];
133 >        double  f0[4], f1[4], hp1, hp2;
134  
135          if (n == 0) {
136                  f[A] = rand3a(xlim[0][i&1],xlim[1][i>>1&1],xlim[2][i>>2]);
# Line 135 | Line 141 | register int  i, n;
141                  n--;
142                  interpolate(f0, i, n);
143                  interpolate(f1, i | 1<<n, n);
144 <                f[A] = (1.0-xarg[n])*f0[A] + xarg[n]*f1[A];
145 <                f[B] = (1.0-xarg[n])*f0[B] + xarg[n]*f1[B];
146 <                f[C] = (1.0-xarg[n])*f0[C] + xarg[n]*f1[C];
147 <                f[D] = hermite(f0[D], f1[D], f0[n], f1[n], xarg[n]);
144 >                hp1 = hpoly1(xarg[n]); hp2 = hpoly2(xarg[n]);
145 >                f[A] = f0[A]*hp1 + f1[A]*hp2;
146 >                f[B] = f0[B]*hp1 + f1[B]*hp2;
147 >                f[C] = f0[C]*hp1 + f1[C]*hp2;
148 >                f[D] = f0[D]*hp1 + f1[D]*hp2 +
149 >                                f0[n]*hpoly3(xarg[n]) + f1[n]*hpoly4(xarg[n]);
150          }
151   }
152  
153  
154 < double
155 < frand(s)                        /* get random number from seed */
156 < register long  s;
154 > static double
155 > frand(                  /* get random number from seed */
156 >        register long  s
157 > )
158   {
159          s = s<<13 ^ s;
160          return(1.0-((s*(s*s*15731+789221)+1376312589)&0x7fffffff)/1073741824.0);
161   }
162  
163  
164 < double
165 < l_hermite()                     /* library call for hermite interpolation */
164 > static double
165 > fnoise3(                        /* compute fractal noise function */
166 >        double  p[3]
167 > )
168   {
169 <        double  t;
159 <        
160 <        t = argument(5);
161 <        return( hermite(argument(1), argument(2),
162 <                        argument(3), argument(4), t) );
163 < }
164 <
165 <
166 < double
167 < fnoise3(p)                      /* compute fractal noise function */
168 < double  p[3];
169 < {
170 <        double  floor();
171 <        long  t[3], v[3], beg[3], s;
169 >        long  t[3], v[3], beg[3];
170          double  fval[8], fc;
171          int  branch;
172 +        register long  s;
173          register int  i, j;
174                                                  /* get starting cube */
175          s = (long)(1.0/EPSILON);
# Line 188 | Line 187 | double  p[3];
187          }
188                                                  /* compute fractal */
189          for ( ; ; ) {
190 <                s >>= 1;
190 >                fc = 0.0;
191 >                for (j = 0; j < 8; j++)
192 >                        fc += fval[j];
193 >                fc *= 0.125;
194 >                if ((s >>= 1) == 0)
195 >                        return(fc);             /* close enough */
196                  branch = 0;
197                  for (i = 0; i < 3; i++) {       /* do center */
198                          v[i] = beg[i] + s;
# Line 196 | Line 200 | double  p[3];
200                                  branch |= 1<<i;
201                          }
202                  }
199                fc = 0.0;
200                for (j = 0; j < 8; j++)
201                        fc += fval[j];
202                fc *= 0.125;
203                if (s < 1)
204                        return(fc);             /* close enough */
203                  fc += s*EPSILON*frand3(v[0],v[1],v[2]);
204                  fval[~branch & 7] = fc;
205                  for (i = 0; i < 3; i++) {       /* do faces */
# Line 218 | Line 216 | double  p[3];
216                          v[i] = beg[i] + s;
217                  }
218                  for (i = 0; i < 3; i++) {       /* do edges */
219 <                        j = (i+1)%3;
219 >                        if ((j = i+1) == 3) j = 0;
220                          if (branch & 1<<j)
221                                  v[j] += s;
222                          else
223                                  v[j] -= s;
224 <                        j = (i+2)%3;
224 >                        if (++j == 3) j = 0;
225                          if (branch & 1<<j)
226                                  v[j] += s;
227                          else
# Line 232 | Line 230 | double  p[3];
230                          fc += fval[branch | 1<<i];
231                          fc = 0.5*fc + s*EPSILON*frand3(v[0],v[1],v[2]);
232                          fval[branch^1<<i] = fc;
233 <                        j = (i+1)%3;
233 >                        if ((j = i+1) == 3) j = 0;
234                          v[j] = beg[j] + s;
235 <                        j = (i+2)%3;
235 >                        if (++j == 3) j = 0;
236                          v[j] = beg[j] + s;
237                  }
238                  for (i = 0; i < 3; i++)         /* new cube */

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines