1 /**************************************************************************
2 * pattern.cpp
3 *
4 * This module implements texturing functions that return a value to be
5 * used in a pigment or normal.
6 *
7 * from Persistence of Vision(tm) Ray Tracer version 3.6.
8 * Copyright 1991-2003 Persistence of Vision Team
9 * Copyright 2003-2004 Persistence of Vision Raytracer Pty. Ltd.
10 *---------------------------------------------------------------------------
11 * NOTICE: This source code file is provided so that users may experiment
12 * with enhancements to POV-Ray and to port the software to platforms other
13 * than those supported by the POV-Ray developers. There are strict rules
14 * regarding how you are permitted to use this file. These rules are contained
15 * in the distribution and derivative versions licenses which should have been
16 * provided with this file.
17 *
18 * These licences may be found online, linked from the end-user license
19 * agreement that is located at http://www.povray.org/povlegal.html
20 *---------------------------------------------------------------------------
21 * This program is based on the popular DKB raytracer version 2.12.
22 * DKBTrace was originally written by David K. Buck.
23 * DKBTrace Ver 2.0-2.12 were written by David K. Buck & Aaron A. Collins.
24 *---------------------------------------------------------------------------
25 *
26 *===========================================================================
27 * This file is part of MegaPOV, a modified and unofficial version of POV-Ray
28 * For more information on MegaPOV visit our website:
29 * http://megapov.inetart.net/
30 *===========================================================================
31 *
32 * $RCSfile: pattern.cpp,v $
33 * $Revision: 1.23 $
34 * $Author: chris $
35 *
36 *****************************************************************************/
37
38 /*
39 * Some texture ideas garnered from SIGGRAPH '85 Volume 19 Number 3,
40 * "An Image Synthesizer" By Ken Perlin.
41 * Further Ideas Garnered from "The RenderMan Companion" (Addison Wesley).
42 */
43
44 #include "frame.h"
45 #include "vector.h"
46 #include "matrices.h"
47 #include "pattern.h"
48 #include "povray.h"
49 #include "texture.h"
50 #include "image.h"
51 #include "txttest.h"
52 #include "colour.h"
53 #include "isosurf.h"
54 #include "pov_util.h"
55 #include "pigment.h"
56 #include "fnpovfpu.h"
57 #include "objects.h"
58 #ifdef PROJECTION_PATCH
59 #include "ray.h"
60 #endif
61 #include <algorithm>
62
63 BEGIN_POV_NAMESPACE
64
65 /*****************************************************************************
66 * Local preprocessor defines
67 ******************************************************************************/
68
69 #define CLIP_DENSITY(r) { if((r) < 0.0) { (r) = 1.0; } else if((r) > 1.0) { (r) = 0.0; } else { (r) = 1.0 - (r); } }
70
71 const int FRACTAL_MAX_EXPONENT = 33;
72
73
74 /*****************************************************************************
75 * Local variables
76 ******************************************************************************/
77
78 bool BinomialCoefficientsInited = false; // GLOBAL VARIABLE
79 int BinomialCoefficients[((FRACTAL_MAX_EXPONENT+1)*(FRACTAL_MAX_EXPONENT+2))/2]; // GLOBAL VARIABLE
80
81
82 /*****************************************************************************
83 * Static functions
84 ******************************************************************************/
85
86 static DBL agate_pattern (VECTOR EPoint, TPATTERN *TPat);
87 static DBL boxed_pattern (VECTOR EPoint);
88 static DBL brick_pattern (VECTOR EPoint, TPATTERN *TPat);
89 static DBL cells_pattern (VECTOR EPoint);
90 static DBL checker_pattern (VECTOR EPoint);
91 static DBL crackle_pattern (VECTOR EPoint, TPATTERN *TPat);
92 static DBL cylindrical_pattern (VECTOR EPoint);
93 static DBL dents_pattern (VECTOR EPoint, TPATTERN *TPat);
94 static DBL density_pattern (VECTOR EPoint, TPATTERN *TPat);
95 static DBL function_pattern (VECTOR EPoint, TPATTERN *TPat); // iso_surface - added
96 static DBL gradient_pattern (VECTOR EPoint, TPATTERN *TPat);
97 static DBL granite_pattern (VECTOR EPoint, TPATTERN *TPat);
98 static DBL hexagon_pattern (VECTOR EPoint);
99 static DBL julia_pattern (VECTOR EPoint, TPATTERN *TPat);
100 static DBL julia3_pattern (VECTOR EPoint, TPATTERN *TPat);
101 static DBL julia4_pattern (VECTOR EPoint, TPATTERN *TPat);
102 static DBL juliax_pattern (VECTOR EPoint, TPATTERN *TPat);
103 static DBL leopard_pattern (VECTOR EPoint);
104 static DBL magnet1m_pattern (VECTOR EPoint, TPATTERN *TPat);
105 static DBL magnet1j_pattern (VECTOR EPoint, TPATTERN *TPat);
106 static DBL magnet2m_pattern (VECTOR EPoint, TPATTERN *TPat);
107 static DBL magnet2j_pattern (VECTOR EPoint, TPATTERN *TPat);
108 static DBL mandel_pattern (VECTOR EPoint, TPATTERN *TPat);
109 static DBL mandel3_pattern (VECTOR EPoint, TPATTERN *TPat);
110 static DBL mandel4_pattern (VECTOR EPoint, TPATTERN *TPat);
111 static DBL mandelx_pattern (VECTOR EPoint, TPATTERN *TPat);
112 static DBL marble_pattern (VECTOR EPoint, TPATTERN *TPat);
113 static DBL object_pattern (VECTOR EPoint, TPATTERN *TPat);
114 static DBL onion_pattern (VECTOR EPoint);
115 static DBL pigment_pattern (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *isect);
116 static DBL planar_pattern (VECTOR EPoint);
117 static DBL quilted_pattern (VECTOR EPoint, TPATTERN *TPat);
118 static DBL radial_pattern (VECTOR EPoint);
119 static DBL ripples_pattern (VECTOR EPoint, TPATTERN *TPat);
120 static DBL slope_pattern (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *Intersection);
121 static DBL spiral1_pattern (VECTOR EPoint, TPATTERN *TPat);
122 static DBL spiral2_pattern (VECTOR EPoint, TPATTERN *TPat);
123 static DBL spherical_pattern (VECTOR EPoint);
124 static DBL waves_pattern (VECTOR EPoint, TPATTERN *TPat);
125 static DBL wood_pattern (VECTOR EPoint, TPATTERN *TPat);
126 static DBL wrinkles_pattern (VECTOR EPoint, TPATTERN *TPat);
127
128 #ifndef REMOVE_NOT_USED_CODE_WARNINGS_PATCH
129 static DBL object(VECTOR EPoint, TPATTERN *TPat);
130 #endif
131
132 #ifdef LISTED_PATTERN_PATCH
133 static DBL ListedPat(VECTOR EPoint, TPATTERN *TPat);
134 #endif
135
136 static DBL fractal_exterior_color(TPATTERN *TPat, int iters, DBL a, DBL b);
137 static DBL fractal_interior_color(TPATTERN *TPat, int iters, DBL a, DBL b, DBL mindist2);
138 static TURB *Search_For_Turb(WARP *Warps);
139 #ifndef REMOVE_NOT_USED_CODE_WARNINGS_PATCH
140 static TURB *Copy_Turb(TURB *Old); // Unused function [AED]
141 #endif
142 static unsigned short readushort(IStream *infile);
143 static unsigned int readuint(IStream *infile);
144 static void InitializeBinomialCoefficients();
145
146 #ifdef ANGLE_OF_INCIDENCE_PATCH
147 static DBL angle_of_incidence (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *Intersection);
148 #endif
149
150 #ifdef PROJECTION_PATCH
151 static DBL projection_pat (VECTOR EPoint, TPATTERN *TPat,INTERSECTION *Inter);
152 #endif
153 /*****************************************************************************
154 *
155 * FUNCTION
156 *
157 * Evaluate_Pattern
158 *
159 * INPUT
160 *
161 * EPoint -- The point in 3d space at which the pattern
162 * is evaluated.
163 * TPat -- Texture pattern struct
164 * Intersection - intersection structure
165 *
166 * OUTPUT
167 *
168 * RETURNS
169 *
170 * DBL result usual 0.0 to 1.0 but may be 2.0 in hexagon
171 *
172 * AUTHOR
173 *
174 * Adapted from Add_Pigment by Chris Young
175 *
176 * DESCRIPTION
177 *
178 * CHANGES
179 *
180 ******************************************************************************/
181
Evaluate_TPat(TPATTERN * TPat,VECTOR EPoint,INTERSECTION * Isection)182 DBL Evaluate_TPat (TPATTERN *TPat, VECTOR EPoint, INTERSECTION *Isection)
183 {
184 DBL value = 0.0;
185
186 /* NK 19 Nov 1999 removed Warp_EPoint call */
187
188 switch(TPat->Type)
189 {
190 case AGATE_PATTERN: value = agate_pattern (EPoint, TPat); break;
191 case BOZO_PATTERN:
192 case SPOTTED_PATTERN:
193 case BUMPS_PATTERN: value = Noise (EPoint, TPat); break;
194 case BRICK_PATTERN: value = brick_pattern (EPoint, TPat); break;
195 case CELLS_PATTERN: value = cells_pattern (EPoint); break;
196 case CHECKER_PATTERN: value = checker_pattern (EPoint); break;
197 case CRACKLE_PATTERN: value = crackle_pattern (EPoint, TPat); break;
198 case GRADIENT_PATTERN: value = gradient_pattern (EPoint, TPat); break;
199 case GRANITE_PATTERN: value = granite_pattern (EPoint, TPat); break;
200 case HEXAGON_PATTERN: value = hexagon_pattern (EPoint); break;
201 case JULIA_PATTERN: value = julia_pattern (EPoint, TPat); break;
202 case JULIA3_PATTERN: value = julia3_pattern (EPoint, TPat); break;
203 case JULIA4_PATTERN: value = julia4_pattern (EPoint, TPat); break;
204 case JULIAX_PATTERN: value = juliax_pattern (EPoint, TPat); break;
205 case LEOPARD_PATTERN: value = leopard_pattern (EPoint); break;
206 case MAGNET1M_PATTERN: value = magnet1m_pattern (EPoint, TPat); break;
207 case MAGNET1J_PATTERN: value = magnet1j_pattern (EPoint, TPat); break;
208 case MAGNET2M_PATTERN: value = magnet2m_pattern (EPoint, TPat); break;
209 case MAGNET2J_PATTERN: value = magnet2j_pattern (EPoint, TPat); break;
210 case MANDEL_PATTERN: value = mandel_pattern (EPoint, TPat); break;
211 case MANDEL3_PATTERN: value = mandel3_pattern (EPoint, TPat); break;
212 case MANDEL4_PATTERN: value = mandel4_pattern (EPoint, TPat); break;
213 case MANDELX_PATTERN: value = mandelx_pattern (EPoint, TPat); break;
214 case MARBLE_PATTERN: value = marble_pattern (EPoint, TPat); break;
215 case ONION_PATTERN: value = onion_pattern (EPoint); break;
216 case RADIAL_PATTERN: value = radial_pattern (EPoint); break;
217 case SPIRAL1_PATTERN: value = spiral1_pattern (EPoint, TPat); break;
218 case SPIRAL2_PATTERN: value = spiral2_pattern (EPoint, TPat); break;
219 case WOOD_PATTERN: value = wood_pattern (EPoint, TPat); break;
220 case WAVES_PATTERN: value = waves_pattern (EPoint, TPat); break;
221 case RIPPLES_PATTERN: value = ripples_pattern (EPoint, TPat); break;
222 case WRINKLES_PATTERN: value = wrinkles_pattern (EPoint, TPat); break;
223 case DENTS_PATTERN: value = dents_pattern (EPoint, TPat); break;
224 case QUILTED_PATTERN: value = quilted_pattern (EPoint, TPat); break;
225 case FUNCTION_PATTERN: value = function_pattern (EPoint, TPat); break;
226 case PLANAR_PATTERN: value = planar_pattern (EPoint); break;
227 case BOXED_PATTERN: value = boxed_pattern (EPoint); break;
228 case SPHERICAL_PATTERN: value = spherical_pattern (EPoint); break;
229 case CYLINDRICAL_PATTERN: value = cylindrical_pattern(EPoint); break;
230 case DENSITY_FILE_PATTERN:value = density_pattern (EPoint, TPat); break;
231 case IMAGE_PATTERN: value = image_pattern (EPoint, TPat); break;
232 case SLOPE_PATTERN: value = slope_pattern (EPoint, TPat, Isection); break;
233 case PIGMENT_PATTERN: value = pigment_pattern (EPoint, TPat, Isection); break;
234 case OBJECT_PATTERN: value = object_pattern (EPoint, TPat); break;
235
236 #ifdef LISTED_PATTERN_PATCH
237 case LISTED_PATTERN: value = ListedPat (EPoint, TPat); break;
238 #endif
239 #ifdef ANGLE_OF_INCIDENCE_PATCH
240 case AOI_PATTERN: value = angle_of_incidence (EPoint, TPat, Isection); break;
241 #endif
242 #ifdef PROJECTION_PATCH
243 case PROJECTION_PATTERN: value = projection_pat (EPoint, TPat,Isection); break;
244 #endif
245 default: Error("Problem in Evaluate_TPat.");
246 }
247
248 if(TPat->Frequency != 0.0)
249 value = fmod(value * TPat->Frequency + TPat->Phase, 1.00001);
250
251 /* allow negative Frequency */
252 if(value < 0.0)
253 value -= floor(value);
254
255 switch(TPat->Wave_Type)
256 {
257 case RAMP_WAVE:
258 break;
259 case SINE_WAVE:
260 value = (1.0 + cycloidal(value)) * 0.5;
261 break;
262 case TRIANGLE_WAVE:
263 value = Triangle_Wave(value);
264 break;
265 case SCALLOP_WAVE:
266 value = fabs(cycloidal(value * 0.5));
267 break;
268 case CUBIC_WAVE:
269 value = Sqr(value) * ((-2.0 * value) + 3.0);
270 break;
271 case POLY_WAVE:
272 value = pow(value, (DBL) TPat->Exponent);
273 break;
274 default:
275 Error("Unknown Wave Type %d.", TPat->Wave_Type);
276 }
277
278 return value;
279 }
280
281
282 /*****************************************************************************
283 *
284 * FUNCTION
285 *
286 * INPUT
287 *
288 * OUTPUT
289 *
290 * RETURNS
291 *
292 * AUTHOR
293 *
294 * DESCRIPTION
295 *
296 * CHANGES
297 *
298 ******************************************************************************/
299
Init_TPat_Fields(TPATTERN * Tpat)300 void Init_TPat_Fields (TPATTERN *Tpat)
301 {
302 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
303 memset(Tpat, 0, sizeof(TPATTERN));
304 #endif
305 Tpat->Type = NO_PATTERN;
306 Tpat->Wave_Type = RAMP_WAVE;
307 Tpat->Flags = NO_FLAGS;
308 Tpat->References = 1;
309 Tpat->Exponent = 1.0;
310 Tpat->Frequency = 1.0;
311 Tpat->Phase = 0.0;
312 Tpat->Warps = NULL;
313 Tpat->Next = NULL;
314 Tpat->Blend_Map = NULL;
315 }
316
317
318 /*****************************************************************************
319 *
320 * FUNCTION
321 *
322 * INPUT
323 *
324 * OUTPUT
325 *
326 * RETURNS
327 *
328 * AUTHOR
329 *
330 * DESCRIPTION
331 *
332 * CHANGES
333 *
334 ******************************************************************************/
335
Copy_TPat_Fields(TPATTERN * New,TPATTERN * Old)336 void Copy_TPat_Fields (TPATTERN *New, TPATTERN *Old)
337 {
338 *New = *Old;
339
340 /* Copy warp chain */
341 New->Warps = Copy_Warps(Old->Warps);
342
343 New->Blend_Map = Copy_Blend_Map(Old->Blend_Map);
344
345 /* Note, cannot copy Old->Next because we don't know what kind of
346 thing this is. It must be copied by Copy_Pigment, Copy_Tnormal etc.
347 */
348
349 /* NK 1998 - added IMAGE_PATTERN */
350 if ((Old->Type == BITMAP_PATTERN) || (Old->Type == IMAGE_PATTERN))
351 {
352 New->Vals.Image = Copy_Image(Old->Vals.Image);
353 }
354
355 if (Old->Type == DENSITY_FILE_PATTERN)
356 {
357 New->Vals.Density_File = Copy_Density_File(Old->Vals.Density_File);
358 }
359
360 if (Old->Type == PIGMENT_PATTERN )
361 {
362 New->Vals.Pigment = Copy_Pigment(Old->Vals.Pigment);
363 }
364
365 if (Old->Type == OBJECT_PATTERN)
366 {
367 if(Old->Vals.Object != NULL)
368 {
369 New->Vals.Object = (OBJECT*)Copy_Object(Old->Vals.Object);
370 }
371 }
372
373 if (Old->Type == CRACKLE_PATTERN)
374 {
375 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
376 New->Vals.Crackle = (CRACKLEPTR)POV_MALLOC(sizeof (CRACKLE), "crackle");
377 *New->Vals.Crackle=*Old->Vals.Crackle;
378 #endif
379 if (Old->Vals.Crackle LESS_MEMORY_PATCH_MO cv != NULL)
380 {
381 New->Vals.Crackle LESS_MEMORY_PATCH_MO cv =(VECTOR*) POV_MALLOC( 125*sizeof(VECTOR), "crackle cache");
382 New->Vals.Crackle LESS_MEMORY_PATCH_MO lastseed = 0x8000000;
383 }
384 }
385
386 if (Old->Type == FUNCTION_PATTERN)
387 {
388 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
389 New->Vals.Function = (FUNCTIONSPTR)POV_MALLOC(sizeof (FUNCTIONS), "functions");
390 *New->Vals.Function=*Old->Vals.Function;
391 #endif
392 if (Old->Vals.Function LESS_MEMORY_PATCH_MO Fn != NULL)
393 {
394 New->Vals.Function LESS_MEMORY_PATCH_MO Fn = (void *)Copy_Function( (FUNCTION_PTR)(Old->Vals.Function LESS_MEMORY_PATCH_MO Fn) );
395 }
396 }
397
398 #ifdef PIGMENT_CAMERA_VIEW_PATCH
399 if (Old->Type == CAMERA_VIEW_PATTERN)
400 {
401 if(Old->Vals.Camera != NULL)
402 {
403 New->Vals.Camera=Copy_Camera( Old->Vals.Camera);
404 }
405 }
406 #endif
407
408 #ifdef PROJECTION_PATCH
409 if (Old->Type == PROJECTION_PATTERN)
410 {
411 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
412 New->Vals.projection = (PROJECTIONPTR)POV_MALLOC(sizeof (PROJECTION), "projection");
413 *New->Vals.projection=*Old->Vals.projection;
414 #endif
415
416 if(Old->Vals.projection LESS_MEMORY_PATCH_MO obj != NULL)
417 {
418 New->Vals.projection LESS_MEMORY_PATCH_MO obj = (OBJECT*)Copy_Object(Old->Vals.projection LESS_MEMORY_PATCH_MO obj);
419 }
420 }
421 #endif
422
423 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
424 #ifdef ANGLE_OF_INCIDENCE_PATCH
425 if (Old->Type == AOI_PATTERN)
426 {
427 New->Vals.Aoi = (AOIPTR)POV_MALLOC(sizeof (AOI), "aoi");
428 *New->Vals.Aoi=*Old->Vals.Aoi;
429 }
430 #endif
431
432 #ifdef NOISE_PIGMENT_PATCH
433 if (Old->Type == NOISE_PIGMENT)
434 {
435 New->Vals.NoisePigment = (NOISEPIGMENTPTR)POV_MALLOC(sizeof (NOISEPIGMENT), "noisepigment");
436 *New->Vals.NoisePigment=*Old->Vals.NoisePigment;
437 }
438 #endif
439
440 if (Old->Type >= JULIA_PATTERN && Old->Type <=MAGNET2J_PATTERN)
441 {
442 New->Vals.Fractal = (FRACTALVALSPTR)POV_MALLOC(sizeof (FRACTALVALS), "fractals");
443 *New->Vals.Fractal=*Old->Vals.Fractal;
444 }
445
446 if (Old->Type ==SLOPE_PATTERN)
447 {
448 New->Vals.Slope = (SLOPEPTR)POV_MALLOC(sizeof (SLOPE), "slope");
449 *New->Vals.Slope=*Old->Vals.Slope;
450 }
451
452 if (Old->Type ==GRADIENT_PATTERN)
453 {
454 New->Vals.Gradient = (VECTOR*)POV_MALLOC(sizeof (VECTOR), "gradient");
455 Assign_Vector(New->Vals.Gradient[0],Old->Vals.Gradient[0]);
456 }
457
458 if (Old->Type ==BRICK_PATTERN)
459 {
460 New->Vals.Brick = (BRICKPTR)POV_MALLOC(sizeof (BRICK), "brick");
461 *New->Vals.Brick=*Old->Vals.Brick;
462 }
463
464 if (Old->Type ==FACETS_PATTERN)
465 {
466 New->Vals.Facets = (FACETSPTR)POV_MALLOC(sizeof (FACETS), "brick");
467 *New->Vals.Facets=*Old->Vals.Facets;
468 }
469
470 if (Old->Type ==QUILTED_PATTERN)
471 {
472 New->Vals.Quilted = (QUILTEDPTR)POV_MALLOC(sizeof (QUILTED), "brick");
473 *New->Vals.Quilted=*Old->Vals.Quilted;
474 }
475 #endif
476 }
477
478
479 /*****************************************************************************
480 *
481 * FUNCTION
482 *
483 * INPUT
484 *
485 * OUTPUT
486 *
487 * RETURNS
488 *
489 * AUTHOR
490 *
491 * DESCRIPTION
492 *
493 * CHANGES
494 *
495 ******************************************************************************/
496
Destroy_TPat_Fields(TPATTERN * Tpat)497 void Destroy_TPat_Fields(TPATTERN *Tpat)
498 {
499 Destroy_Warps(Tpat->Warps);
500 Destroy_Blend_Map(Tpat->Blend_Map);
501 /* Note, cannot destroy Tpat->Next nor pattern itself because we don't
502 know what kind of thing this is. It must be destroied by Destroy_Pigment, etc.
503 */
504
505 if ((Tpat->Type == BITMAP_PATTERN) || (Tpat->Type == IMAGE_PATTERN))
506 {
507 Destroy_Image(Tpat->Vals.Image);
508 }
509
510 if (Tpat->Type == DENSITY_FILE_PATTERN)
511 {
512 Destroy_Density_File(Tpat->Vals.Density_File);
513 }
514
515 if (Tpat->Type == OBJECT_PATTERN)
516 {
517 if(Tpat->Vals.Object != NULL)
518 {
519 Destroy_Object((OBJECT *)Tpat->Vals.Object);
520 }
521 }
522
523 if (Tpat->Type == CRACKLE_PATTERN)
524 {
525 if (Tpat->Vals.Crackle LESS_MEMORY_PATCH_MO cv != NULL)
526 {
527 POV_FREE( Tpat->Vals.Crackle LESS_MEMORY_PATCH_MO cv );
528 }
529 }
530
531 if (Tpat->Type == PIGMENT_PATTERN)
532 {
533 if (Tpat->Vals.Pigment != NULL)
534 {
535 Destroy_Pigment( Tpat->Vals.Pigment );
536 }
537 }
538
539 if (Tpat->Type == FUNCTION_PATTERN)
540 {
541 if (Tpat->Vals.Function LESS_MEMORY_PATCH_MO Fn != NULL)
542 {
543 Destroy_Function( (FUNCTION_PTR)(Tpat->Vals.Function LESS_MEMORY_PATCH_MO Fn) );
544 }
545 }
546
547 #ifdef PIGMENT_CAMERA_VIEW_PATCH
548 if (Tpat->Type == CAMERA_VIEW_PATTERN)
549 {
550 if(Tpat->Vals.Camera != NULL)
551 {
552 Destroy_Camera((CAMERA *)Tpat->Vals.Camera );
553 }
554 }
555 #endif
556
557 #ifdef PROJECTION_PATCH
558 if (Tpat->Type == PROJECTION_PATTERN)
559 {
560 if(Tpat->Vals.projection LESS_MEMORY_PATCH_MO obj != NULL)
561 {
562 Destroy_Object((OBJECT *)Tpat->Vals.projection LESS_MEMORY_PATCH_MO obj);
563 }
564 }
565 #endif
566
567 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
568 #ifdef PROJECTION_PATCH
569 if (Tpat->Type == PROJECTION_PATTERN)
570 {
571 if(Tpat->Vals.projection != NULL)
572 {
573 POV_FREE( Tpat->Vals.projection );
574 }
575 }
576 #endif
577
578 #ifdef ANGLE_OF_INCIDENCE_PATCH
579 if (Tpat->Type == AOI_PATTERN)
580 {
581 if(Tpat->Vals.Aoi != NULL)
582 {
583 POV_FREE( Tpat->Vals.Aoi );
584 }
585 }
586 #endif
587
588 #ifdef NOISE_PIGMENT_PATCH
589 if (Tpat->Type == NOISE_PIGMENT)
590 {
591 if(Tpat->Vals.NoisePigment != NULL)
592 {
593 POV_FREE( Tpat->Vals.NoisePigment );
594 }
595 }
596 #endif
597
598
599 if (Tpat->Type >= JULIA_PATTERN && Tpat->Type <=MAGNET2J_PATTERN)
600 {
601 if(Tpat->Vals.Fractal != NULL)
602 {
603 POV_FREE( Tpat->Vals.Fractal );
604 }
605 }
606
607 if (Tpat->Type == CRACKLE_PATTERN )
608 {
609 if(Tpat->Vals.Crackle != NULL)
610 {
611 POV_FREE( Tpat->Vals.Crackle );
612 }
613 }
614
615 if (Tpat->Type ==SLOPE_PATTERN)
616 {
617 if(Tpat->Vals.Slope != NULL)
618 {
619 POV_FREE( Tpat->Vals.Slope );
620 }
621 }
622
623 if (Tpat->Type ==GRADIENT_PATTERN)
624 {
625 if(Tpat->Vals.Gradient != NULL)
626 {
627 POV_FREE( Tpat->Vals.Gradient );
628 }
629 }
630
631 if (Tpat->Type ==BRICK_PATTERN)
632 {
633 if(Tpat->Vals.Brick != NULL)
634 {
635 POV_FREE( Tpat->Vals.Brick );
636 }
637 }
638
639 if (Tpat->Type ==FACETS_PATTERN)
640 {
641 if(Tpat->Vals.Facets != NULL)
642 {
643 POV_FREE( Tpat->Vals.Facets );
644 }
645 }
646
647 if (Tpat->Type ==QUILTED_PATTERN)
648 {
649 if(Tpat->Vals.Quilted != NULL)
650 {
651 POV_FREE( Tpat->Vals.Quilted );
652 }
653 }
654
655 if (Tpat->Type ==FUNCTION_PATTERN)
656 {
657 if(Tpat->Vals.Function != NULL)
658 {
659 POV_FREE( Tpat->Vals.Function );
660 }
661 }
662
663 #endif
664 }
665
666
667 /*****************************************************************************
668 *
669 * FUNCTION
670 *
671 * INPUT
672 *
673 * OUTPUT
674 *
675 * RETURNS
676 *
677 * AUTHOR
678 *
679 * DESCRIPTION
680 *
681 * CHANGES
682 *
683 ******************************************************************************/
684
Create_Turb()685 TURB *Create_Turb()
686 {
687 TURB *New;
688
689 New = (TURB *)POV_MALLOC(sizeof(TURB),"turbulence struct");
690
691 Make_Vector(New->Turbulence, 0.0, 0.0, 0.0);
692
693 New->Octaves = 6;
694 New->Omega = 0.5;
695 New->Lambda = 2.0;
696
697 return(New);
698 }
699
700
701 /*****************************************************************************
702 *
703 * FUNCTION
704 *
705 * INPUT
706 *
707 * OUTPUT
708 *
709 * RETURNS
710 *
711 * AUTHOR
712 *
713 * DESCRIPTION
714 *
715 * CHANGES
716 *
717 ******************************************************************************/
718
719 #if 0 /* Unused function [AED] */
720 static TURB *Copy_Turb(TURB *Old)
721 {
722 TURB *New;
723
724 if (Old != NULL)
725 {
726 New = Create_Turb();
727
728 *New = *Old;
729 }
730 else
731 {
732 New=NULL;
733 }
734
735 return(New);
736 }
737 #endif
738
739
740 /*****************************************************************************
741 *
742 * FUNCTION
743 *
744 * Translate_Tpattern
745 *
746 * INPUT
747 *
748 * OUTPUT
749 *
750 * RETURNS
751 *
752 * AUTHOR
753 *
754 * POV-Ray Team
755 *
756 * DESCRIPTION
757 *
758 * CHANGES
759 *
760 ******************************************************************************/
761
Translate_Tpattern(TPATTERN * Tpattern,VECTOR Vector)762 void Translate_Tpattern(TPATTERN *Tpattern,VECTOR Vector)
763 {
764 TRANSFORM Trans;
765
766 if (Tpattern != NULL)
767 {
768 Compute_Translation_Transform (&Trans, Vector);
769
770 Transform_Tpattern (Tpattern, &Trans);
771 }
772 }
773
774
775 /*****************************************************************************
776 *
777 * FUNCTION
778 *
779 * Rotate_Tpattern
780 *
781 * INPUT
782 *
783 * OUTPUT
784 *
785 * RETURNS
786 *
787 * AUTHOR
788 *
789 * POV-Ray Team
790 *
791 * DESCRIPTION
792 *
793 * CHANGES
794 *
795 ******************************************************************************/
796
Rotate_Tpattern(TPATTERN * Tpattern,VECTOR Vector)797 void Rotate_Tpattern(TPATTERN *Tpattern,VECTOR Vector)
798 {
799 TRANSFORM Trans;
800
801 if (Tpattern != NULL)
802 {
803 Compute_Rotation_Transform (&Trans, Vector);
804
805 Transform_Tpattern (Tpattern, &Trans);
806 }
807 }
808
809
810 /*****************************************************************************
811 *
812 * FUNCTION
813 *
814 * Scale_Tpattern
815 *
816 * INPUT
817 *
818 * OUTPUT
819 *
820 * RETURNS
821 *
822 * AUTHOR
823 *
824 * POV-Ray Team
825 *
826 * DESCRIPTION
827 *
828 * CHANGES
829 *
830 ******************************************************************************/
831
Scale_Tpattern(TPATTERN * Tpattern,VECTOR Vector)832 void Scale_Tpattern(TPATTERN *Tpattern,VECTOR Vector)
833 {
834 TRANSFORM Trans;
835
836 if (Tpattern != NULL)
837 {
838 Compute_Scaling_Transform (&Trans, Vector);
839
840 Transform_Tpattern (Tpattern, &Trans);
841 }
842 }
843
844
845 /*****************************************************************************
846 *
847 * FUNCTION
848 *
849 * Transform_Tpattern
850 *
851 * INPUT
852 *
853 * OUTPUT
854 *
855 * RETURNS
856 *
857 * AUTHOR
858 *
859 * POV-Ray Team
860 *
861 * DESCRIPTION
862 *
863 * CHANGES
864 *
865 ******************************************************************************/
866
Transform_Tpattern(TPATTERN * Tpattern,TRANSFORM * Trans)867 void Transform_Tpattern(TPATTERN *Tpattern,TRANSFORM *Trans)
868 {
869 WARP *Temp;
870
871 if (Tpattern != NULL)
872 {
873 if (Tpattern->Warps == NULL)
874 {
875 Tpattern->Warps = Create_Warp(TRANSFORM_WARP);
876 }
877 else
878 {
879 if (Tpattern->Warps->Warp_Type != TRANSFORM_WARP)
880 {
881 Temp = Tpattern->Warps;
882
883 Tpattern->Warps = Create_Warp(TRANSFORM_WARP);
884
885 Tpattern->Warps->Next_Warp = Temp;
886 if(Tpattern->Warps->Next_Warp != NULL)
887 Tpattern->Warps->Next_Warp->Prev_Warp = Tpattern->Warps;
888 }
889 }
890
891 Compose_Transforms (&( ((TRANS *)(Tpattern->Warps))->Trans), Trans);
892 }
893 }
894
895
896 /*****************************************************************************
897 *
898 * FUNCTION
899 *
900 * INPUT
901 *
902 * OUTPUT
903 *
904 * RETURNS
905 *
906 * AUTHOR
907 *
908 * DESCRIPTION
909 *
910 * CHANGES
911 *
912 ******************************************************************************/
913
Search_Blend_Map(DBL value,BLEND_MAP * Blend_Map,BLEND_MAP_ENTRY ** Prev,BLEND_MAP_ENTRY ** Cur)914 void Search_Blend_Map (DBL value,BLEND_MAP *Blend_Map,BLEND_MAP_ENTRY **Prev,BLEND_MAP_ENTRY **Cur)
915 {
916 BLEND_MAP_ENTRY *P, *C;
917 int Max_Ent=Blend_Map->Number_Of_Entries-1;
918
919 /* if greater than last, use last. */
920
921 if (value >= Blend_Map->Blend_Map_Entries[Max_Ent].value)
922 {
923 P = C = &(Blend_Map->Blend_Map_Entries[Max_Ent]);
924 }
925 else
926 {
927 P = C = &(Blend_Map->Blend_Map_Entries[0]);
928
929 while (value > C->value)
930 {
931 P = C++;
932 }
933 }
934
935 if (value == C->value)
936 {
937 P = C;
938 }
939
940 *Prev = P;
941 *Cur = C;
942 }
943
944
945 /*****************************************************************************
946 *
947 * FUNCTION
948 *
949 * INPUT
950 *
951 * OUTPUT
952 *
953 * RETURNS
954 *
955 * AUTHOR
956 *
957 * DESCRIPTION
958 *
959 * CHANGES
960 *
961 ******************************************************************************/
962
Search_For_Turb(WARP * Warps)963 static TURB *Search_For_Turb(WARP *Warps)
964 {
965 WARP* Temp=Warps;
966
967 if (Temp!=NULL)
968 {
969 while (Temp->Next_Warp != NULL)
970 {
971 Temp=Temp->Next_Warp;
972 }
973
974 if (Temp->Warp_Type != CLASSIC_TURB_WARP)
975 {
976 Temp=NULL;
977 }
978 }
979
980 return ((TURB *)Temp);
981 }
982
983
984 /*****************************************************************************
985 *
986 * FUNCTION
987 *
988 * agate_pattern
989 *
990 * INPUT
991 *
992 * EPoint -- The point in 3d space at which the pattern is evaluated.
993 * TPat -- Texture pattern struct
994 *
995 * OUTPUT
996 *
997 * RETURNS
998 *
999 * DBL value in the range 0.0 to 1.0
1000 *
1001 * AUTHOR
1002 *
1003 * POV-Ray Team
1004 *
1005 * DESCRIPTION
1006 *
1007 * CHANGES
1008 *
1009 * Oct 1994 : adapted from agate pigment by [CY]
1010 *
1011 ******************************************************************************/
1012
agate_pattern(VECTOR EPoint,TPATTERN * TPat)1013 static DBL agate_pattern (VECTOR EPoint, TPATTERN *TPat)
1014 {
1015 register DBL noise, turb_val;
1016 TURB* Turb;
1017
1018 Turb=Search_For_Turb(TPat->Warps);
1019
1020 turb_val = TPat->Vals.Agate_Turb_Scale * Turbulence(EPoint,Turb,TPat);
1021
1022 noise = 0.5 * (cycloidal(1.3 * turb_val + 1.1 * EPoint[Z]) + 1.0);
1023
1024 if (noise < 0.0)
1025 {
1026 noise = 0.0;
1027 }
1028 else
1029 {
1030 noise = min(1.0, noise);
1031 noise = pow(noise, 0.77);
1032 }
1033
1034 return(noise);
1035 }
1036
1037
1038 /*****************************************************************************
1039 *
1040 * FUNCTION
1041 *
1042 * boxed_pattern
1043 *
1044 * INPUT
1045 *
1046 * EPoint -- The point in 3d space at which the pattern
1047 * is evaluated.
1048 *
1049 * OUTPUT
1050 *
1051 * RETURNS
1052 *
1053 * DBL value in the range 0.0 to 1.0
1054 *
1055 * AUTHOR
1056 *
1057 * -
1058 *
1059 * DESCRIPTION
1060 *
1061 * -
1062 *
1063 * CHANGES
1064 *
1065 * -
1066 *
1067 ******************************************************************************/
1068
boxed_pattern(VECTOR EPoint)1069 static DBL boxed_pattern (VECTOR EPoint)
1070 {
1071 register DBL value;
1072
1073 value = max(fabs(EPoint[X]), max(fabs(EPoint[Y]), fabs(EPoint[Z])));
1074 CLIP_DENSITY(value);
1075
1076 return(value);
1077 }
1078
1079
1080 /*****************************************************************************
1081 *
1082 * FUNCTION
1083 *
1084 * brick_pattern
1085 *
1086 * INPUT
1087 *
1088 * EPoint -- The point in 3d space at which the pattern
1089 * is evaluated.
1090 * TPat -- Texture pattern struct
1091 *
1092 * OUTPUT
1093 *
1094 * RETURNS
1095 *
1096 * DBL value exactly 0.0 or 1.0
1097 *
1098 * AUTHOR
1099 *
1100 * Dan Farmer
1101 *
1102 * DESCRIPTION
1103 *
1104 * CHANGES
1105 *
1106 * Oct 1994 : adapted from pigment by [CY]
1107 *
1108 ******************************************************************************/
1109
brick_pattern(VECTOR EPoint,TPATTERN * TPat)1110 static DBL brick_pattern (VECTOR EPoint, TPATTERN *TPat)
1111 {
1112 int ibrickx, ibricky, ibrickz;
1113 DBL brickheight, brickwidth, brickdepth;
1114 DBL brickmortar, mortarheight, mortarwidth, mortardepth;
1115 DBL brickx, bricky, brickz;
1116 DBL x, y, z, fudgit;
1117
1118 fudgit=Small_Tolerance+TPat->Vals.Brick LESS_MEMORY_PATCH_MO Mortar;
1119
1120 x = EPoint[X]+fudgit;
1121 y = EPoint[Y]+fudgit;
1122 z = EPoint[Z]+fudgit;
1123
1124 brickwidth = TPat->Vals.Brick LESS_MEMORY_PATCH_MO Size[X];
1125 brickheight = TPat->Vals.Brick LESS_MEMORY_PATCH_MO Size[Y];
1126 brickdepth = TPat->Vals.Brick LESS_MEMORY_PATCH_MO Size[Z];
1127 brickmortar = (DBL)TPat->Vals.Brick LESS_MEMORY_PATCH_MO Mortar;
1128
1129 mortarwidth = brickmortar / brickwidth;
1130 mortarheight = brickmortar / brickheight;
1131 mortardepth = brickmortar / brickdepth;
1132
1133 /* 1) Check mortar layers in the X-Z plane (ie: top view) */
1134
1135 bricky = y / brickheight;
1136 ibricky = (int) bricky;
1137 bricky -= (DBL) ibricky;
1138
1139 if (bricky < 0.0)
1140 {
1141 bricky += 1.0;
1142 }
1143
1144 if (bricky <= mortarheight)
1145 {
1146 return(0.0);
1147 }
1148
1149 bricky = (y / brickheight) * 0.5;
1150 ibricky = (int) bricky;
1151 bricky -= (DBL) ibricky;
1152
1153 if (bricky < 0.0)
1154 {
1155 bricky += 1.0;
1156 }
1157
1158
1159 /* 2) Check ODD mortar layers in the Y-Z plane (ends) */
1160
1161 brickx = (x / brickwidth);
1162 ibrickx = (int) brickx;
1163 brickx -= (DBL) ibrickx;
1164
1165 if (brickx < 0.0)
1166 {
1167 brickx += 1.0;
1168 }
1169
1170 if ((brickx <= mortarwidth) && (bricky <= 0.5))
1171 {
1172 return(0.0);
1173 }
1174
1175 /* 3) Check EVEN mortar layers in the Y-Z plane (ends) */
1176
1177 brickx = (x / brickwidth) + 0.5;
1178 ibrickx = (int) brickx;
1179 brickx -= (DBL) ibrickx;
1180
1181 if (brickx < 0.0)
1182 {
1183 brickx += 1.0;
1184 }
1185
1186 if ((brickx <= mortarwidth) && (bricky > 0.5))
1187 {
1188 return(0.0);
1189 }
1190
1191 /* 4) Check ODD mortar layers in the Y-X plane (facing) */
1192
1193 brickz = (z / brickdepth);
1194 ibrickz = (int) brickz;
1195 brickz -= (DBL) ibrickz;
1196
1197 if (brickz < 0.0)
1198 {
1199 brickz += 1.0;
1200 }
1201
1202 if ((brickz <= mortardepth) && (bricky > 0.5))
1203 {
1204 return(0.0);
1205 }
1206
1207 /* 5) Check EVEN mortar layers in the X-Y plane (facing) */
1208
1209 brickz = (z / brickdepth) + 0.5;
1210 ibrickz = (int) brickz;
1211 brickz -= (DBL) ibrickz;
1212
1213 if (brickz < 0.0)
1214 {
1215 brickz += 1.0;
1216 }
1217
1218 if ((brickz <= mortardepth) && (bricky <= 0.5))
1219 {
1220 return(0.0);
1221 }
1222
1223 /* If we've gotten this far, color me brick. */
1224
1225 return(1.0);
1226 }
1227
1228
1229 /*****************************************************************************
1230 *
1231 * FUNCTION
1232 *
1233 * cells_pattern
1234 *
1235 * INPUT
1236 *
1237 * EPoint -- The point in 3d space at which the pattern
1238 * is evaluated.
1239 * OUTPUT
1240 *
1241 * RETURNS
1242 *
1243 * DBL value in the range 0.0 to 1.0
1244 *
1245 * AUTHOR
1246 *
1247 * John VanSickle
1248 *
1249 * DESCRIPTION
1250 *
1251 * "cells":
1252 *
1253 * New colour function by John VanSickle,
1254 * vansickl@erols.com
1255 *
1256 * Assigns a pseudorandom value to each unit cube. The value for the cube in
1257 * which the evaluted point lies is returned.
1258 *
1259 * All "cells" specific source code and examples are in the public domain.
1260 *
1261 * CHANGES
1262 *
1263 * -
1264 *
1265 ******************************************************************************/
1266
cells_pattern(VECTOR EPoint)1267 static DBL cells_pattern (VECTOR EPoint)
1268 {
1269 int temp;
1270 DBL tf;
1271
1272 temp = POV_GET_OLD_RAND(); /* save current seed */
1273
1274
1275 /* select a random value based on the cube from which this came. */
1276
1277 /* floor the values, instead of just truncating - this eliminates duplicated cells
1278 around the axes */
1279
1280 POV_SRAND(Hash3d((int)floor(EPoint[X]+Small_Tolerance),
1281 (int)floor(EPoint[Y]+Small_Tolerance),
1282 (int)floor(EPoint[Z]+Small_Tolerance)));
1283
1284 tf = FRAND();
1285
1286 POV_SRAND(temp); /* restore */
1287
1288 return min(tf, 1.0);
1289 }
1290
1291
1292 /*****************************************************************************
1293 *
1294 * FUNCTION
1295 *
1296 * checker_pattern
1297 *
1298 * INPUT
1299 *
1300 * EPoint -- The point in 3d space at which the pattern
1301 * is evaluated.
1302 *
1303 * OUTPUT
1304 *
1305 * RETURNS
1306 *
1307 * DBL value exactly 0.0 or 1.0
1308 *
1309 * AUTHOR
1310 *
1311 * POV-Team
1312 *
1313 * DESCRIPTION
1314 *
1315 * CHANGES
1316 * Oct 1994 : adapted from pigment by [CY]
1317 *
1318 ******************************************************************************/
1319
checker_pattern(VECTOR EPoint)1320 static DBL checker_pattern (VECTOR EPoint)
1321 {
1322 int value;
1323
1324 value = (int)(floor(EPoint[X]+Small_Tolerance) +
1325 floor(EPoint[Y]+Small_Tolerance) +
1326 floor(EPoint[Z]+Small_Tolerance));
1327
1328 if (value & 1)
1329 {
1330 return (1.0);
1331 }
1332 else
1333 {
1334 return (0.0);
1335 }
1336 }
1337
1338
1339 /*****************************************************************************
1340 *
1341 * FUNCTION
1342 *
1343 * crackle_pattern
1344 *
1345 * INPUT
1346 *
1347 * EPoint -- The point in 3d space at which the pattern
1348 * is evaluated.
1349 * OUTPUT
1350 *
1351 * RETURNS
1352 *
1353 * DBL value in the range 0.0 to 1.0
1354 *
1355 * AUTHOR
1356 *
1357 * Jim McElhiney
1358 *
1359 * DESCRIPTION
1360 *
1361 * "crackle":
1362 *
1363 * New colour function by Jim McElhiney,
1364 * CompuServe 71201,1326, aka mcelhiney@acm.org
1365 *
1366 * Large scale, without turbulence, makes a pretty good stone wall.
1367 * Small scale, without turbulence, makes a pretty good crackle ceramic glaze.
1368 * Highly turbulent (with moderate displacement) makes a good marble, solving
1369 * the problem of apparent parallel layers in Perlin's method.
1370 * 2 octaves of full-displacement turbulence make a great "drizzled paint"
1371 * pattern, like a 1950's counter top.
1372 * Rule of thumb: put a single colour transition near 0 in your colour map.
1373 *
1374 * Mathematically, the set crackle(p)=0 is a 3D Voronoi diagram of a field of
1375 * semirandom points, and crackle(p)>0 is distance from set along shortest path.
1376 * (A Voronoi diagram is the locus of points equidistant from their 2 nearest
1377 * neighbours from a set of disjoint points, like the membranes in suds are
1378 * to the centres of the bubbles).
1379 *
1380 * All "crackle" specific source code and examples are in the public domain.
1381 *
1382 * CHANGES
1383 * Oct 1994 : adapted from pigment by [CY]
1384 * Other changes: enhanced by Ron Parker, Integer math by Nathan Kopp
1385 *
1386 ******************************************************************************/
1387 static int IntPickInCube(int tvx, int tvy, int tvz, VECTOR p1);
1388
crackle_pattern(VECTOR EPoint,TPATTERN * TPat)1389 static DBL crackle_pattern (VECTOR EPoint, TPATTERN *TPat )
1390 {
1391 int i;
1392 int thisseed;
1393 DBL sum, minsum, minsum2, minsum3, tf;
1394 VECTOR minvec;
1395 VECTOR tv, dv, t1;
1396 int addx,addy,addz;
1397
1398 VECTOR flo;
1399 int cvc;
1400 #ifdef REMOVE_NOT_USED_CODE_WARNINGS_PATCH
1401 static int vals[3];
1402 #else
1403 static int vali=0, vals[3];
1404 #endif
1405 static int valid[125];
1406
1407 DBL Metric;
1408 DBL Offset;
1409
1410 int UseSquare;
1411 int UseUnity;
1412 int flox,floy,floz;
1413 /*int seed,temp;*/
1414
1415 Metric = TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Metric[X];
1416 Offset = TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Offset;
1417
1418 UseSquare = ( Metric == 2);
1419 UseUnity = ( Metric == 1);
1420
1421 Assign_Vector(tv,EPoint);
1422
1423 /*
1424 * Check to see if the input point is in the same unit cube as the last
1425 * call to this function, to use cache of cubelets for speed.
1426 */
1427
1428 thisseed = PickInCube(tv, t1);
1429
1430 if (true)//(thisseed != TPat->Vals.Crackle LESS_MEMORY_PATCH_MO lastseed)
1431 {
1432 /*
1433 * No, not same unit cube. Calculate the random points for this new
1434 * cube and its 80 neighbours which differ in any axis by 1 or 2.
1435 * Why distance of 2? If there is 1 point in each cube, located
1436 * randomly, it is possible for the closest random point to be in the
1437 * cube 2 over, or the one two over and one up. It is NOT possible
1438 * for it to be two over and two up. Picture a 3x3x3 cube with 9 more
1439 * cubes glued onto each face.
1440 */
1441
1442
1443 flo[X] = floor(tv[X] - EPSILON);
1444 flo[Y] = floor(tv[Y] - EPSILON);
1445 flo[Z] = floor(tv[Z] - EPSILON);
1446
1447 Assign_Vector( (TPat->Vals.Crackle LESS_MEMORY_PATCH_MO lastcenter), flo );
1448
1449 /* Now store a points for this cube and each of the 80 neighbour cubes. */
1450
1451 vals[0]=25*(2+(-2))+5*(2+(-1))+2+(-1);
1452 vals[1]=25*(2+(-2))+5*(2+(-1))+2+(0);
1453 vals[2]=25*(2+(-2))+5*(2+(-1))+2+(1);
1454
1455 flox = (int)flo[X];
1456 floy = (int)flo[Y];
1457 floz = (int)flo[Z];
1458
1459 for (addx = -2; addx <= 2; addx++)
1460 {
1461 for (addy = -2; addy <= 2; addy++)
1462 {
1463 for (addz = -2; addz <= 2; addz++)
1464 {
1465 /* For each cubelet in a 5x5 cube. */
1466 cvc = 25*(2+addx)+5*(2+addy)+2+addz;
1467
1468 if ((abs(addx)==2)+(abs(addy)==2)+(abs(addz)==2) <= 1)
1469 {
1470 /* Yes, it's within a 3d knight move away. */
1471
1472 #define INLINE_PICK_IN_CUBE 0
1473 #if INLINE_PICK_IN_CUBE
1474 /* do our own PickInCube and use as much integer math as possible */
1475 seed = Hash3d((flox+addx)&0xFFF,(floy+addy)&0xFFF,(floz+addz)&0xFFF);
1476 temp = POV_GET_OLD_RAND(); /* save current seed */
1477 POV_SRAND(seed);
1478 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][X] = flox+addx + FRAND();
1479 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][Y] = floy+addy + FRAND();
1480 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][Z] = floz+addz + FRAND();
1481 POV_SRAND(temp); /* restore */
1482 #else
1483 IntPickInCube(flox+addx,floy+addy,floz+addz, t1);
1484
1485 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][X] = t1[X];
1486 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][Y] = t1[Y];
1487 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[cvc][Z] = t1[Z];
1488 #endif
1489 valid[cvc]=1;
1490 }
1491 else
1492 {
1493 valid[cvc]=0;
1494 }
1495 }
1496 }
1497 }
1498
1499 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO lastseed = thisseed;
1500 }
1501
1502 cvc=125;
1503 /*
1504 * Find the 2 points with the 2 shortest distances from the input point.
1505 * Loop invariant: minsum is shortest dist, minsum2 is 2nd shortest
1506 */
1507
1508 /* Set up the loop so the invariant is true: minsum <= minsum2 */
1509
1510 VSub(dv, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[0]], tv);
1511 if ( UseSquare )
1512 {
1513 minsum = VSumSqr(dv);
1514 if ( Offset ) minsum += Offset*Offset;
1515 }
1516 else if ( UseUnity )
1517 {
1518 minsum = fabs(dv[X]) + fabs(dv[Y]) + fabs(dv[Z]);
1519 if ( Offset ) minsum += Offset;
1520 }
1521 else
1522 {
1523 minsum = pow( fabs( dv[X] ), Metric ) +
1524 pow( fabs( dv[Y] ), Metric ) +
1525 pow( fabs( dv[Z] ), Metric );
1526 if ( Offset ) minsum += pow( Offset, Metric );
1527 }
1528 Assign_Vector( minvec, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[0]] );
1529 VSub(dv, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[1]], tv);
1530 if ( UseSquare )
1531 {
1532 minsum2 = VSumSqr(dv);
1533 if ( Offset ) minsum2 += Offset*Offset;
1534 }
1535 else if ( UseUnity )
1536 {
1537 minsum2 = fabs(dv[X]) + fabs(dv[Y]) + fabs(dv[Z]);
1538 if ( Offset ) minsum2 += Offset;
1539 }
1540 else
1541 {
1542 minsum2 = pow( fabs( dv[X] ), Metric ) +
1543 pow( fabs( dv[Y] ), Metric ) +
1544 pow( fabs( dv[Z] ), Metric );
1545 if ( Offset ) minsum2 += pow( Offset, Metric );
1546 }
1547 VSub(dv, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[2]], tv);
1548 if ( UseSquare ) {
1549 minsum3 = VSumSqr(dv);
1550 if ( Offset ) minsum3 += Offset*Offset;
1551 }
1552 else if ( UseUnity )
1553 {
1554 minsum3 = fabs(dv[X]) + fabs(dv[Y]) + fabs(dv[Z]);
1555 if ( Offset ) minsum3 += Offset;
1556 }
1557 else
1558 {
1559 minsum3 = pow( fabs( dv[X] ), Metric ) +
1560 pow( fabs( dv[Y] ), Metric ) +
1561 pow( fabs( dv[Z] ), Metric );
1562 if ( Offset ) minsum3 += pow( Offset, Metric );
1563 }
1564
1565 if (minsum2 < minsum)
1566 {
1567 tf = minsum; minsum = minsum2; minsum2 = tf;
1568 Assign_Vector( minvec, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[1]] );
1569 }
1570 if (minsum3 < minsum)
1571 {
1572 tf = minsum; minsum = minsum3; minsum3 = tf;
1573 Assign_Vector( minvec, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[vals[2]] );
1574 }
1575 if ( minsum3 < minsum2 )
1576 {
1577 tf = minsum2; minsum2=minsum3; minsum3= tf;
1578 }
1579
1580 /* Loop for the 81 cubelets to find closest and 2nd closest. */
1581
1582 for (i = vals[2]+1; i < cvc; i++) if (valid[i])
1583 {
1584 VSub(dv, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[i], tv);
1585
1586 if ( UseSquare )
1587 {
1588 sum = VSumSqr(dv);
1589 if ( Offset ) sum += Offset*Offset;
1590 }
1591 else if ( UseUnity )
1592 {
1593 sum = fabs(dv[X]) + fabs(dv[Y]) + fabs(dv[Z]);
1594 if ( Offset ) sum += Offset;
1595 }
1596 else
1597 {
1598 sum = pow( fabs( dv[X] ), Metric ) +
1599 pow( fabs( dv[Y] ), Metric ) +
1600 pow( fabs( dv[Z] ), Metric );
1601 if ( Offset ) sum += pow( Offset, Metric );
1602 }
1603
1604 if (sum < minsum)
1605 {
1606 minsum3 = minsum2;
1607 minsum2 = minsum;
1608 minsum = sum;
1609 Assign_Vector( minvec, TPat->Vals.Crackle LESS_MEMORY_PATCH_MO cv[i] );
1610 }
1611 else if (sum < minsum2)
1612 {
1613 minsum3 = minsum2;
1614 minsum2 = sum;
1615 }
1616 else if ( sum < minsum3 )
1617 {
1618 minsum3 = sum;
1619 }
1620 }
1621
1622 if ( TPat->Vals.Crackle LESS_MEMORY_PATCH_MO IsSolid )
1623 {
1624 tf = Noise( minvec, TPat );
1625 }
1626 else if (UseSquare)
1627 {
1628 tf = TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[X]*sqrt(minsum) +
1629 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Y]*sqrt(minsum2) +
1630 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Z]*sqrt(minsum3);
1631 }
1632 else if ( UseUnity )
1633 {
1634 tf = TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[X]*minsum +
1635 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Y]*minsum2 +
1636 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Z]*minsum3;
1637 }
1638 else
1639 {
1640 tf = TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[X]*pow(minsum, 1.0/Metric) +
1641 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Y]*pow(minsum2, 1.0/Metric) +
1642 TPat->Vals.Crackle LESS_MEMORY_PATCH_MO Form[Z]*pow(minsum3, 1.0/Metric);
1643 }
1644
1645 return max(min(tf, 1.), 0.);
1646 }
1647
1648
1649 /*****************************************************************************
1650 *
1651 * FUNCTION
1652 *
1653 * cylindrical_pattern
1654 *
1655 * INPUT
1656 *
1657 * EPoint -- The point in 3d space at which the pattern
1658 * is evaluated.
1659 *
1660 * OUTPUT
1661 *
1662 * RETURNS
1663 *
1664 * DBL value in the range 0.0 to 1.0
1665 *
1666 * AUTHOR
1667 *
1668 * -
1669 *
1670 * DESCRIPTION
1671 *
1672 * -
1673 *
1674 * CHANGES
1675 *
1676 * -
1677 *
1678 ******************************************************************************/
1679
cylindrical_pattern(VECTOR EPoint)1680 static DBL cylindrical_pattern (VECTOR EPoint)
1681 {
1682 register DBL value;
1683
1684 value = sqrt(Sqr(EPoint[X]) + Sqr(EPoint[Z]));
1685 CLIP_DENSITY(value);
1686
1687 return(value);
1688 }
1689
1690
1691 /*****************************************************************************
1692 *
1693 * FUNCTION
1694 *
1695 * density_pattern
1696 *
1697 * INPUT
1698 *
1699 * OUTPUT
1700 *
1701 * RETURNS
1702 *
1703 * AUTHOR
1704 *
1705 * Dieter Bayer
1706 *
1707 * DESCRIPTION
1708 *
1709 * CHANGES
1710 *
1711 * Dec 1996 : Creation.
1712 *
1713 ******************************************************************************/
1714
intp3(float t,float fa,float fb,float fc,float fd)1715 inline float intp3(float t, float fa, float fb, float fc, float fd)
1716 {
1717 float b,d,e,f;
1718
1719 b = (fc - fa) * 0.5;
1720 d = (fd - fb) * 0.5;
1721 e = 2.0 * (fb - fc) + b + d;
1722 f = -3.0 * (fb - fc) - 2.0 * b - d;
1723
1724 return ((e * t + f) * t + b) * t + fb;
1725 }
1726
intp3_2(float t,float fa,float fb,float fc,float fd)1727 inline float intp3_2(float t, float fa, float fb, float fc, float fd)
1728 {
1729 float b,e,f;
1730
1731 e = fd - fc - fa + fb;
1732 f = fa - fb - e;
1733 b = fc - fa;
1734
1735 return ((e * t + f) * t + b) * t + fb;
1736 }
1737
1738 #define zmax(i,imax) (((i)<0)?(imax-1):((i) % (imax)))
1739
density_pattern(VECTOR EPoint,TPATTERN * TPat)1740 static DBL density_pattern(VECTOR EPoint, TPATTERN *TPat)
1741 {
1742 int x, y, z;
1743 int x1, y1, z1;
1744 int x2, y2, z2;
1745 DBL Ex, Ey, Ez;
1746 DBL xx, yy, zz;
1747 DBL xi, yi;
1748 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
1749 DBL f111=0.0, f112=0.0, f121=0.0, f122=0.0, f211=0.0, f212=0.0, f221=0.0, f222=0.0;
1750 #else
1751 DBL f111, f112, f121, f122, f211, f212, f221, f222;
1752 #endif
1753 float intpd2[4][4];
1754 DBL density = 0.0;
1755 DENSITY_FILE_DATA *Data;
1756 int k0, k1, k2, k3, i,j,ii,jj;
1757
1758 Ex=EPoint[X];
1759 Ey=EPoint[Y];
1760 Ez=EPoint[Z];
1761
1762 if((TPat->Vals.Density_File != NULL) && ((Data = TPat->Vals.Density_File->Data) != NULL))
1763 {
1764 /* if(Data->Cyclic == true)
1765 {
1766 Ex -= floor(Ex);
1767 Ey -= floor(Ey);
1768 Ez -= floor(Ez);
1769 }
1770 */
1771 if((Ex >= 0.0) && (Ex < 1.0) && (Ey >= 0.0) && (Ey < 1.0) && (Ez >= 0.0) && (Ez < 1.0))
1772 {
1773 switch (TPat->Vals.Density_File->Interpolation % 10)
1774 {
1775 case NO_INTERPOLATION:
1776 x = (int)(Ex * (DBL)Data->Sx);
1777 y = (int)(Ey * (DBL)Data->Sy);
1778 z = (int)(Ez * (DBL)Data->Sz);
1779
1780 if ((x < 0) || (x >= Data->Sx) || (y < 0) || (y >= Data->Sy) || (z < 0) || (z >= Data->Sz))
1781 density = 0.0;
1782 else
1783 {
1784 if(Data->Type == 4)
1785 density = (DBL)Data->Density32[z * Data->Sy * Data->Sx + y * Data->Sx + x] / (DBL)UINT_MAX;
1786 else if(Data->Type==2)
1787 density = (DBL)Data->Density16[z * Data->Sy * Data->Sx + y * Data->Sx + x] / (DBL)USHRT_MAX;
1788 else if(Data->Type == 1)
1789 density = (DBL)Data->Density8[z * Data->Sy * Data->Sx + y * Data->Sx + x] / (DBL)UCHAR_MAX;
1790 }
1791 break;
1792 case TRILINEAR_INTERPOLATION:
1793 xx = Ex * (DBL)(Data->Sx );
1794 yy = Ey * (DBL)(Data->Sy );
1795 zz = Ez * (DBL)(Data->Sz );
1796
1797 x1 = (int)xx;
1798 y1 = (int)yy;
1799 z1 = (int)zz;
1800
1801 x2 = (x1 + 1) % Data->Sx;
1802 y2 = (y1 + 1) % Data->Sy;
1803 z2 = (z1 + 1) % Data->Sz;
1804
1805 xx -= floor(xx);
1806 yy -= floor(yy);
1807 zz -= floor(zz);
1808
1809 xi = 1.0 - xx;
1810 yi = 1.0 - yy;
1811
1812 if(Data->Type == 4)
1813 {
1814 f111 = (DBL)Data->Density32[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)UINT_MAX;
1815 f112 = (DBL)Data->Density32[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)UINT_MAX;
1816 f121 = (DBL)Data->Density32[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)UINT_MAX;
1817 f122 = (DBL)Data->Density32[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)UINT_MAX;
1818 f211 = (DBL)Data->Density32[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)UINT_MAX;
1819 f212 = (DBL)Data->Density32[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)UINT_MAX;
1820 f221 = (DBL)Data->Density32[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)UINT_MAX;
1821 f222 = (DBL)Data->Density32[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)UINT_MAX;
1822 }
1823 else if(Data->Type == 2)
1824 {
1825 f111 = (DBL)Data->Density16[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)USHRT_MAX;
1826 f112 = (DBL)Data->Density16[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)USHRT_MAX;
1827 f121 = (DBL)Data->Density16[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)USHRT_MAX;
1828 f122 = (DBL)Data->Density16[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)USHRT_MAX;
1829 f211 = (DBL)Data->Density16[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)USHRT_MAX;
1830 f212 = (DBL)Data->Density16[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)USHRT_MAX;
1831 f221 = (DBL)Data->Density16[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)USHRT_MAX;
1832 f222 = (DBL)Data->Density16[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)USHRT_MAX;
1833 }
1834 else if(Data->Type == 1)
1835 {
1836 f111 = (DBL)Data->Density8[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)UCHAR_MAX;
1837 f112 = (DBL)Data->Density8[z1 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)UCHAR_MAX;
1838 f121 = (DBL)Data->Density8[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)UCHAR_MAX;
1839 f122 = (DBL)Data->Density8[z1 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)UCHAR_MAX;
1840 f211 = (DBL)Data->Density8[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x1] / (DBL)UCHAR_MAX;
1841 f212 = (DBL)Data->Density8[z2 * Data->Sy * Data->Sx + y1 * Data->Sx + x2] / (DBL)UCHAR_MAX;
1842 f221 = (DBL)Data->Density8[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x1] / (DBL)UCHAR_MAX;
1843 f222 = (DBL)Data->Density8[z2 * Data->Sy * Data->Sx + y2 * Data->Sx + x2] / (DBL)UCHAR_MAX;
1844 }
1845
1846 density = ((f111 * xi + f112 * xx) * yi + (f121 * xi + f122 * xx) * yy) * (1.0 - zz) +
1847 ((f211 * xi + f212 * xx) * yi + (f221 * xi + f222 * xx) * yy) * zz;
1848 break;
1849 case TRICUBIC_INTERPOLATION:
1850 default:
1851 xx = Ex * (DBL)(Data->Sx);
1852 yy = Ey * (DBL)(Data->Sy);
1853 zz = Ez * (DBL)(Data->Sz);
1854
1855 x1 = (int)xx;
1856 y1 = (int)yy;
1857 z1 = (int)zz;
1858
1859 xx -= floor(xx);
1860 yy -= floor(yy);
1861 zz -= floor(zz);
1862
1863 k0 = zmax(-1+z1, Data->Sz );
1864 k1 = zmax( z1, Data->Sz );
1865 k2 = zmax( 1+z1, Data->Sz );
1866 k3 = zmax( 2+z1, Data->Sz );
1867
1868 if(Data->Type == 4)
1869 {
1870 for(i = 0; i < 4; i++)
1871 {
1872 ii = zmax(i + x1 - 1, Data->Sx);
1873 for(j = 0; j < 4; j++)
1874 {
1875 jj = zmax(j + y1 - 1, Data->Sy);
1876 intpd2[i][j] = intp3(zz,
1877 Data->Density32[k0 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UINT_MAX,
1878 Data->Density32[k1 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UINT_MAX,
1879 Data->Density32[k2 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UINT_MAX,
1880 Data->Density32[k3 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UINT_MAX);
1881 }
1882 }
1883 }
1884 else if(Data->Type == 2)
1885 {
1886 for(i = 0; i < 4; i++)
1887 {
1888 ii = zmax(i + x1 - 1, Data->Sx);
1889 for(j = 0; j < 4; j++)
1890 {
1891 jj = zmax(j + y1 - 1, Data->Sy);
1892 intpd2[i][j] = intp3(zz,
1893 Data->Density16[k0 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)USHRT_MAX,
1894 Data->Density16[k1 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)USHRT_MAX,
1895 Data->Density16[k2 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)USHRT_MAX,
1896 Data->Density16[k3 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)USHRT_MAX);
1897 }
1898 }
1899 }
1900 else if(Data->Type == 1)
1901 {
1902 for(i = 0; i < 4; i++)
1903 {
1904 ii = zmax(i + x1 - 1, Data->Sx);
1905 for(j = 0; j < 4; j++)
1906 {
1907 jj = zmax(j + y1 - 1, Data->Sy);
1908 intpd2[i][j] = intp3(zz,
1909 Data->Density8[k0 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UCHAR_MAX,
1910 Data->Density8[k1 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UCHAR_MAX,
1911 Data->Density8[k2 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UCHAR_MAX,
1912 Data->Density8[k3 * Data->Sy * Data->Sx + jj * Data->Sx + ii] / (DBL)UCHAR_MAX);
1913 }
1914 }
1915 }
1916
1917 for(i = 0; i < 4; i++)
1918 intpd2[0][i] = intp3(yy, intpd2[i][0], intpd2[i][1], intpd2[i][2], intpd2[i][3]);
1919
1920 density = intp3(xx, intpd2[0][0], intpd2[0][1], intpd2[0][2], intpd2[0][3]);
1921 break;
1922 }
1923 }
1924 else
1925 density = 0.0;
1926 }
1927
1928 return density;
1929 }
1930
1931
1932 /*****************************************************************************
1933 *
1934 * FUNCTION
1935 *
1936 * dents_pattern
1937 *
1938 * INPUT
1939 *
1940 * EPoint -- The point in 3d space at which the pattern
1941 * is evaluated.
1942 *
1943 * OUTPUT
1944 *
1945 * RETURNS
1946 *
1947 * DBL value in the range 0.0 to 1.0
1948 *
1949 * AUTHOR
1950 *
1951 * POV-Ray Team
1952 *
1953 * DESCRIPTION : Note this pattern is only used for pigments and textures.
1954 * Normals have a specialized pattern for this.
1955 *
1956 * CHANGES
1957 * Nov 1994 : adapted from normal by [CY]
1958 *
1959 ******************************************************************************/
1960
dents_pattern(VECTOR EPoint,TPATTERN * TPat)1961 static DBL dents_pattern (VECTOR EPoint, TPATTERN *TPat)
1962 {
1963 DBL noise;
1964
1965 noise = Noise (EPoint, TPat);
1966
1967 return(noise * noise * noise);
1968 }
1969
1970
1971 /*****************************************************************************
1972 *
1973 * FUNCTION
1974 *
1975 * function_pattern
1976 *
1977 * INPUT
1978 *
1979 * EPoint -- The point in 3d space at which the pattern
1980 * is evaluated.
1981 *
1982 * OUTPUT
1983 *
1984 * RETURNS
1985 *
1986 * DBL value in the range 0.0 to 1.0
1987 *
1988 * AUTHOR
1989 *
1990 * POV-Ray Team
1991 *
1992 * DESCRIPTION
1993 *
1994 * CHANGES
1995 *
1996 ******************************************************************************/
1997
function_pattern(VECTOR EPoint,TPATTERN * TPat)1998 static DBL function_pattern (VECTOR EPoint, TPATTERN *TPat)
1999 {
2000 DBL value;
2001 FPUContext *oldcontext;
2002
2003 if(TPat->Vals.Function LESS_MEMORY_PATCH_MO Data == NULL)
2004 TPat->Vals.Function LESS_MEMORY_PATCH_MO Data = POVFPU_NewContext();
2005
2006 oldcontext = POVFPU_SwitchContext((FPUContext *)(TPat->Vals.Function LESS_MEMORY_PATCH_MO Data));
2007
2008 POVFPU_SetLocal(X, EPoint[X]);
2009 POVFPU_SetLocal(Y, EPoint[Y]);
2010 POVFPU_SetLocal(Z, EPoint[Z]);
2011
2012 value = POVFPU_Run(*((FUNCTION_PTR)(TPat->Vals.Function LESS_MEMORY_PATCH_MO Fn)));
2013
2014 (void)POVFPU_SwitchContext(oldcontext);
2015
2016 return ((value > 1.0) ? fmod(value, 1.0) : value);
2017 }
2018
2019
2020 /*****************************************************************************
2021 *
2022 * FUNCTION
2023 *
2024 * gradient_pattern
2025 *
2026 * INPUT
2027 *
2028 * EPoint -- The point in 3d space at which the pattern
2029 * is evaluated.
2030 *
2031 * OUTPUT
2032 *
2033 * RETURNS
2034 *
2035 * DBL value in the range 0.0 to 1.0
2036 *
2037 * AUTHOR
2038 *
2039 * POV-Ray Team
2040 *
2041 * DESCRIPTION
2042 *
2043 * Gradient Pattern - gradient based on the fractional values of
2044 * x, y or z, based on whether or not the given directional vector is
2045 * a 1.0 or a 0.0.
2046 * The basic concept of this is from DBW Render, but Dave Wecker's
2047 * only supports simple Y axis gradients.
2048 *
2049 * CHANGES
2050 *
2051 * Oct 1994 : adapted from pigment by [CY]
2052 *
2053 ******************************************************************************/
2054
gradient_pattern(VECTOR EPoint,TPATTERN * TPat)2055 static DBL gradient_pattern (VECTOR EPoint, TPATTERN *TPat)
2056 {
2057 DBL Result;
2058 #ifdef LESS_MEMORY_IN_PATTERNS_PATCH
2059 VDot( Result, EPoint, *TPat->Vals.Gradient );
2060 #else
2061 VDot( Result, EPoint, TPat->Vals.Gradient );
2062 #endif
2063 /* Mod to keep within [0.0,1.0] range */
2064 return ((Result > 1.0) ? fmod(Result, 1.0) : Result);
2065 }
2066
2067
2068 /*****************************************************************************
2069 *
2070 * FUNCTION
2071 *
2072 * granite_pattern
2073 *
2074 * INPUT
2075 *
2076 * EPoint -- The point in 3d space at which the pattern
2077 * is evaluated.
2078 *
2079 * OUTPUT
2080 *
2081 * RETURNS
2082 *
2083 * DBL value in the range 0.0 to 1.0
2084 *
2085 * AUTHOR
2086 *
2087 * POV-Ray Team
2088 *
2089 * DESCRIPTION
2090 *
2091 * Granite - kind of a union of the "spotted" and the "dented" textures,
2092 * using a 1/f fractal noise function for color values. Typically used
2093 * with small scaling values. Should work with colour maps for pink granite.
2094 *
2095 * CHANGES
2096 *
2097 * Oct 1994 : adapted from pigment by [CY]
2098 *
2099 ******************************************************************************/
2100
granite_pattern(VECTOR EPoint,TPATTERN * TPat)2101 static DBL granite_pattern (VECTOR EPoint, TPATTERN *TPat)
2102 {
2103 register int i;
2104 register DBL temp, noise = 0.0, freq = 1.0;
2105 VECTOR tv1,tv2;
2106
2107 VScale(tv1,EPoint,4.0);
2108
2109 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2110 int noise_generator = ( TPat != NULL ? (TPat->Flags & NOISE_FLAGS) >> 4 : 0 ) ;
2111 #else
2112 int noise_generator;
2113 if (TPat != NULL)
2114 noise_generator = (TPat->Flags & NOISE_FLAGS) >> 4;
2115 #endif
2116
2117 if (!noise_generator)
2118 noise_generator=opts.Noise_Generator;
2119
2120 for (i = 0; i < 6 ; freq *= 2.0, i++)
2121 {
2122 VScale(tv2,tv1,freq);
2123 if(noise_generator==1)
2124 {
2125 temp = 0.5 - Noise (tv2, TPat);
2126 temp = fabs(temp);
2127 }
2128 else
2129 {
2130 temp = 1.0 - 2.0 * Noise (tv2, TPat);
2131 temp = fabs(temp);
2132 if (temp>0.5) temp=0.5;
2133 }
2134
2135
2136
2137 noise += temp / freq;
2138 }
2139
2140 return(noise);
2141 }
2142
2143
2144 /*****************************************************************************
2145 *
2146 * FUNCTION
2147 *
2148 * hexagon_pattern
2149 *
2150 * INPUT
2151 *
2152 * EPoint -- The point in 3d space at which the pattern
2153 * is evaluated.
2154 *
2155 * OUTPUT
2156 *
2157 * RETURNS
2158 *
2159 * DBL value exactly 0.0, 1.0 or 2.0
2160 *
2161 * AUTHOR
2162 *
2163 * Ernest MacDougal Campbell III
2164 *
2165 * DESCRIPTION
2166 *
2167 * TriHex pattern -- Ernest MacDougal Campbell III (EMC3) 11/23/92
2168 *
2169 * Creates a hexagon pattern in the XZ plane.
2170 *
2171 * This algorithm is hard to explain. First it scales the point to make
2172 * a few of the later calculations easier, then maps some points to be
2173 * closer to the Origin. A small area in the first quadrant is subdivided
2174 * into a 6 x 6 grid. The position of the point mapped into that grid
2175 * determines its color. For some points, just the grid location is enough,
2176 * but for others, we have to calculate which half of the block it's in
2177 * (this is where the atan2() function comes in handy).
2178 *
2179 * CHANGES
2180 *
2181 * Nov 1992 : Creation.
2182 * Oct 1994 : adapted from pigment by [CY]
2183 *
2184 ******************************************************************************/
2185
2186 const DBL xfactor = 0.5; /* each triangle is split in half for the grid */
2187 const DBL zfactor = 0.866025404; /* sqrt(3)/2 -- Height of an equilateral triangle */
2188
hexagon_pattern(VECTOR EPoint)2189 static DBL hexagon_pattern (VECTOR EPoint)
2190 {
2191 int xm, zm;
2192 int brkindx;
2193 DBL xs, zs, xl, zl, value = 0.0;
2194 DBL x=EPoint[X];
2195 DBL z=EPoint[Z];
2196
2197
2198 /* Keep all numbers positive. Also, if z is negative, map it in such a
2199 * way as to avoid mirroring across the x-axis. The value 5.196152424
2200 * is (sqrt(3)/2) * 6 (because the grid is 6 blocks high)
2201 */
2202
2203 x = fabs(x);
2204
2205 /* Avoid mirroring across x-axis. */
2206
2207 z = z < 0.0 ? 5.196152424 - fabs(z) : z;
2208
2209 /* Scale point to make calcs easier. */
2210
2211 xs = x/xfactor;
2212 zs = z/zfactor;
2213
2214 /* Map points into the 6 x 6 grid where the basic formula works. */
2215
2216 xs -= floor(xs/6.0) * 6.0;
2217 zs -= floor(zs/6.0) * 6.0;
2218
2219 /* Get a block in the 6 x 6 grid. */
2220
2221 xm = (int) FLOOR(xs) % 6;
2222 zm = (int) FLOOR(zs) % 6;
2223
2224 switch (xm)
2225 {
2226 /* These are easy cases: Color depends only on xm and zm. */
2227
2228 case 0:
2229 case 5:
2230
2231 switch (zm)
2232 {
2233 case 0:
2234 case 5: value = 0; break;
2235
2236 case 1:
2237 case 2: value = 1; break;
2238
2239 case 3:
2240 case 4: value = 2; break;
2241 }
2242
2243 break;
2244
2245 case 2:
2246 case 3:
2247
2248 switch (zm)
2249 {
2250 case 0:
2251 case 1: value = 2; break;
2252
2253 case 2:
2254 case 3: value = 0; break;
2255
2256 case 4:
2257 case 5: value = 1; break;
2258 }
2259
2260 break;
2261
2262 /* These cases are harder. These blocks are divided diagonally
2263 * by the angled edges of the hexagons. Some slope positive, and
2264 * others negative. We flip the x value of the negatively sloped
2265 * pieces. Then we check to see if the point in question falls
2266 * in the upper or lower half of the block. That info, plus the
2267 * z status of the block determines the color.
2268 */
2269
2270 case 1:
2271 case 4:
2272
2273 /* Map the point into the block at the origin. */
2274
2275 xl = xs-xm;
2276 zl = zs-zm;
2277
2278 /* These blocks have negative slopes so we flip it horizontally. */
2279
2280 if (((xm + zm) % 2) == 1)
2281 {
2282 xl = 1.0 - xl;
2283 }
2284
2285 /* Avoid a divide-by-zero error. */
2286
2287 if (xl == 0.0)
2288 {
2289 xl = 0.0001;
2290 }
2291
2292 /* Is the angle less-than or greater-than 45 degrees? */
2293
2294 brkindx = (zl / xl) < 1.0;
2295
2296 /* was...
2297 * brkindx = (atan2(zl,xl) < (45 * M_PI_180));
2298 * ...but because of the mapping, it's easier and cheaper,
2299 * CPU-wise, to just use a good ol' slope.
2300 */
2301
2302 switch (brkindx)
2303 {
2304 case true:
2305
2306 switch (zm)
2307 {
2308 case 0:
2309 case 3: value = 0; break;
2310
2311 case 2:
2312 case 5: value = 1; break;
2313
2314 case 1:
2315 case 4: value = 2; break;
2316 }
2317
2318 break;
2319
2320 case false:
2321
2322 switch (zm)
2323 {
2324 case 0:
2325 case 3: value = 2; break;
2326
2327 case 2:
2328 case 5: value = 0; break;
2329
2330 case 1:
2331 case 4: value = 1; break;
2332 }
2333
2334 break;
2335 }
2336 }
2337
2338 value = fmod(value, 3.0);
2339
2340 return(value);
2341 }
2342
2343
2344 /*****************************************************************************
2345 *
2346 * FUNCTION
2347 *
2348 * julia_pattern
2349 *
2350 * INPUT
2351 *
2352 * EPoint -- The point in 3d space at which the pattern
2353 * is evaluated.
2354 *
2355 * OUTPUT
2356 *
2357 * RETURNS
2358 *
2359 * DBL value in the range 0.0 to 1.0
2360 *
2361 * AUTHOR
2362 *
2363 * Nieminen Juha
2364 *
2365 * DESCRIPTION
2366 *
2367 * -
2368 *
2369 * CHANGES
2370 *
2371 * -
2372 *
2373 ******************************************************************************/
2374
julia_pattern(VECTOR EPoint,TPATTERN * TPat)2375 static DBL julia_pattern (VECTOR EPoint, TPATTERN *TPat)
2376 {
2377 int it_max, col;
2378 DBL a, b, cf
2379 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2380 =0
2381 #endif
2382 , a2, b2, dist2, mindist2,
2383 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V];
2384
2385 a = EPoint[X]; a2 = Sqr(a);
2386 b = EPoint[Y]; b2 = Sqr(b);
2387 mindist2 = a2+b2;
2388
2389 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2390
2391 for (col = 0; col < it_max; col++)
2392 {
2393 b = 2.0 * a * b + ci;
2394 a = a2 - b2 + cr;
2395
2396 a2 = Sqr(a);
2397 b2 = Sqr(b);
2398 dist2 = a2+b2;
2399
2400 if(dist2 < mindist2) mindist2 = dist2;
2401 if(dist2 > 4.0)
2402 {
2403 cf = fractal_exterior_color(TPat, col, a, b);
2404 break;
2405 }
2406 }
2407
2408 if(col == it_max)
2409 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2410
2411 return(cf);
2412 }
2413
2414
2415 /*****************************************************************************
2416 *
2417 * FUNCTION
2418 *
2419 * julia3_pattern
2420 *
2421 * INPUT
2422 *
2423 * EPoint -- The point in 3d space at which the pattern
2424 * is evaluated.
2425 *
2426 * OUTPUT
2427 *
2428 * RETURNS
2429 *
2430 * DBL value in the range 0.0 to 1.0
2431 *
2432 * AUTHOR
2433 *
2434 * Nieminen Juha
2435 *
2436 * DESCRIPTION
2437 *
2438 * -
2439 *
2440 * CHANGES
2441 *
2442 * -
2443 *
2444 ******************************************************************************/
2445
julia3_pattern(VECTOR EPoint,TPATTERN * TPat)2446 static DBL julia3_pattern (VECTOR EPoint, TPATTERN *TPat)
2447 {
2448 int it_max, col;
2449 DBL a, b, cf
2450 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2451 =0
2452 #endif
2453 , a2, b2, dist2, mindist2,
2454 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V];
2455
2456 a = EPoint[X]; a2 = Sqr(a);
2457 b = EPoint[Y]; b2 = Sqr(b);
2458 mindist2 = a2+b2;
2459
2460 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2461
2462 for (col = 0; col < it_max; col++)
2463 {
2464 b = 3.0*a2*b - b2*b + ci;
2465 a = a2*a - 3.0*a*b2 + cr;
2466
2467 a2 = Sqr(a);
2468 b2 = Sqr(b);
2469 dist2 = a2+b2;
2470
2471 if(dist2 < mindist2) mindist2 = dist2;
2472 if(dist2 > 4.0)
2473 {
2474 cf = fractal_exterior_color(TPat, col, a, b);
2475 break;
2476 }
2477 }
2478
2479 if(col == it_max)
2480 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2481
2482 return(cf);
2483 }
2484
2485
2486 /*****************************************************************************
2487 *
2488 * FUNCTION
2489 *
2490 * julia4_pattern
2491 *
2492 * INPUT
2493 *
2494 * EPoint -- The point in 3d space at which the pattern
2495 * is evaluated.
2496 *
2497 * OUTPUT
2498 *
2499 * RETURNS
2500 *
2501 * DBL value in the range 0.0 to 1.0
2502 *
2503 * AUTHOR
2504 *
2505 * Nieminen Juha
2506 *
2507 * DESCRIPTION
2508 *
2509 * -
2510 *
2511 * CHANGES
2512 *
2513 * -
2514 *
2515 ******************************************************************************/
2516
julia4_pattern(VECTOR EPoint,TPATTERN * TPat)2517 static DBL julia4_pattern (VECTOR EPoint, TPATTERN *TPat)
2518 {
2519 int it_max, col;
2520 DBL a, b, cf
2521 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2522 =0
2523 #endif
2524 , a2, b2, dist2, mindist2,
2525 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V];
2526
2527 a = EPoint[X]; a2 = Sqr(a);
2528 b = EPoint[Y]; b2 = Sqr(b);
2529 mindist2 = a2+b2;
2530
2531 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2532
2533 for (col = 0; col < it_max; col++)
2534 {
2535 b = 4.0 * (a2*a*b - a*b2*b) + ci;
2536 a = a2*a2 - 6.0*a2*b2 + b2*b2 + cr;
2537
2538 a2 = Sqr(a);
2539 b2 = Sqr(b);
2540 dist2 = a2+b2;
2541
2542 if(dist2 < mindist2) mindist2 = dist2;
2543 if(dist2 > 4.0)
2544 {
2545 cf = fractal_exterior_color(TPat, col, a, b);
2546 break;
2547 }
2548 }
2549
2550 if(col == it_max)
2551 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2552
2553 return(cf);
2554 }
2555
2556
2557 /*****************************************************************************
2558 *
2559 * FUNCTION
2560 *
2561 * juliax_pattern
2562 *
2563 * INPUT
2564 *
2565 * EPoint -- The point in 3d space at which the pattern
2566 * is evaluated.
2567 *
2568 * OUTPUT
2569 *
2570 * RETURNS
2571 *
2572 * DBL value in the range 0.0 to 1.0
2573 *
2574 * AUTHOR
2575 *
2576 * Nieminen Juha
2577 *
2578 * DESCRIPTION
2579 *
2580 * -
2581 *
2582 * CHANGES
2583 *
2584 * -
2585 *
2586 ******************************************************************************/
2587
juliax_pattern(VECTOR EPoint,TPATTERN * TPat)2588 static DBL juliax_pattern (VECTOR EPoint, TPATTERN *TPat)
2589 {
2590 int it_max, col, exponent;
2591 DBL a, b, cf=0, x, y, dist2, mindist2,
2592 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V];
2593 int* binomial_coeff;
2594
2595 if(BinomialCoefficientsInited == false)
2596 InitializeBinomialCoefficients();
2597
2598 a = x = EPoint[X];
2599 b = y = EPoint[Y];
2600 mindist2 = a*a+b*b;
2601
2602 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2603 exponent = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Exponent;
2604
2605 binomial_coeff = &BinomialCoefficients[(exponent+1)*exponent/2];
2606
2607 for (col = 0; col < it_max; col++)
2608 {
2609 // Calculate (a+bi)^exponent
2610 DBL new_a = pow(a, exponent);
2611 for(int k=2; k<=exponent; k+=2)
2612 {
2613 new_a += binomial_coeff[k]*pow(a, exponent-k)*pow(b, k);
2614 }
2615 DBL new_b = 0;
2616 for(int l=1; l<=exponent; l+=2)
2617 {
2618 new_b += binomial_coeff[l]*pow(a, exponent-l)*pow(b, l);
2619 }
2620
2621 a = new_a + cr;
2622 b = new_b + ci;
2623
2624 dist2 = a*a+b*b;
2625
2626 if(dist2 < mindist2) mindist2 = dist2;
2627 if(dist2 > 4.0)
2628 {
2629 cf = fractal_exterior_color(TPat, col, a, b);
2630 break;
2631 }
2632 }
2633
2634 if(col == it_max)
2635 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2636
2637 return(cf);
2638 }
2639
2640
2641 /*****************************************************************************
2642 *
2643 * FUNCTION
2644 *
2645 * leopard_pattern
2646 *
2647 * INPUT
2648 *
2649 * EPoint -- The point in 3d space at which the pattern
2650 * is evaluated.
2651 *
2652 * OUTPUT
2653 *
2654 * RETURNS
2655 *
2656 * DBL value in the range 0.0 to 1.0
2657 *
2658 * AUTHOR
2659 *
2660 * Scott Taylor
2661 *
2662 * DESCRIPTION
2663 *
2664 * CHANGES
2665 *
2666 * Jul 1991 : Creation.
2667 * Oct 1994 : adapted from pigment by [CY]
2668 *
2669 ******************************************************************************/
2670
leopard_pattern(VECTOR EPoint)2671 static DBL leopard_pattern (VECTOR EPoint)
2672 {
2673 register DBL value, temp1, temp2, temp3;
2674
2675 /* This form didn't work with Zortech 386 compiler */
2676 /* value = Sqr((sin(x)+sin(y)+sin(z))/3); */
2677 /* So we break it down. */
2678
2679 temp1 = sin(EPoint[X]);
2680 temp2 = sin(EPoint[Y]);
2681 temp3 = sin(EPoint[Z]);
2682
2683 value = Sqr((temp1 + temp2 + temp3) / 3.0);
2684
2685 return(value);
2686 }
2687
2688
2689 /*****************************************************************************
2690 *
2691 * FUNCTION
2692 *
2693 * magnet1m_pattern
2694 *
2695 * INPUT
2696 *
2697 * EPoint -- The point in 3d space at which the pattern
2698 * is evaluated.
2699 *
2700 * OUTPUT
2701 *
2702 * RETURNS
2703 *
2704 * DBL value in the range 0.0 to 1.0
2705 *
2706 * AUTHOR
2707 *
2708 * Nieminen Juha
2709 *
2710 * DESCRIPTION
2711 *
2712 * -
2713 *
2714 * CHANGES
2715 *
2716 * -
2717 *
2718 ******************************************************************************/
2719
magnet1m_pattern(VECTOR EPoint,TPATTERN * TPat)2720 static DBL magnet1m_pattern (VECTOR EPoint, TPATTERN *TPat)
2721 {
2722 int it_max, col;
2723 DBL a, b, cf
2724 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2725 =0
2726 #endif
2727 , a2, b2, x, y, tmp, tmp1r, tmp1i, tmp2r, tmp2i, dist2, mindist2;
2728
2729 x = EPoint[X];
2730 y = EPoint[Y];
2731 a = a2 = 0;
2732 b = b2 = 0;
2733 mindist2 = 10000;
2734
2735 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2736
2737 for (col = 0; col < it_max; col++)
2738 {
2739 tmp1r = a2-b2 + x-1;
2740 tmp1i = 2*a*b + y;
2741 tmp2r = 2*a + x-2;
2742 tmp2i = 2*b + y;
2743 tmp = tmp2r*tmp2r + tmp2i*tmp2i;
2744 a = (tmp1r*tmp2r + tmp1i*tmp2i) / tmp;
2745 b = (tmp1i*tmp2r - tmp1r*tmp2i) / tmp;
2746 b2 = b*b;
2747 b = 2*a*b;
2748 a = a*a-b2;
2749
2750 a2 = Sqr(a);
2751 b2 = Sqr(b);
2752 dist2 = a2+b2;
2753
2754 if(dist2 < mindist2) mindist2 = dist2;
2755 tmp1r = a-1;
2756 if(dist2 > 10000.0 || tmp1r*tmp1r+b2 < 1/10000.0)
2757 {
2758 cf = fractal_exterior_color(TPat, col, a, b);
2759 break;
2760 }
2761 }
2762
2763 if(col == it_max)
2764 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2765
2766 return(cf);
2767 }
2768
2769
2770 /*****************************************************************************
2771 *
2772 * FUNCTION
2773 *
2774 * magnet1j_pattern
2775 *
2776 * INPUT
2777 *
2778 * EPoint -- The point in 3d space at which the pattern
2779 * is evaluated.
2780 *
2781 * OUTPUT
2782 *
2783 * RETURNS
2784 *
2785 * DBL value in the range 0.0 to 1.0
2786 *
2787 * AUTHOR
2788 *
2789 * Nieminen Juha
2790 *
2791 * DESCRIPTION
2792 *
2793 * -
2794 *
2795 * CHANGES
2796 *
2797 * -
2798 *
2799 ******************************************************************************/
2800
magnet1j_pattern(VECTOR EPoint,TPATTERN * TPat)2801 static DBL magnet1j_pattern (VECTOR EPoint, TPATTERN *TPat)
2802 {
2803 int it_max, col;
2804 DBL a, b, cf
2805 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2806 =0
2807 #endif
2808 , a2, b2, tmp, tmp1r, tmp1i, tmp2r, tmp2i, dist2, mindist2,
2809 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V];
2810
2811 a = EPoint[X]; a2 = Sqr(a);
2812 b = EPoint[Y]; b2 = Sqr(b);
2813 mindist2 = a2+b2;
2814
2815 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2816
2817 for (col = 0; col < it_max; col++)
2818 {
2819 tmp1r = a2-b2 + cr-1;
2820 tmp1i = 2*a*b + ci;
2821 tmp2r = 2*a + cr-2;
2822 tmp2i = 2*b + ci;
2823 tmp = tmp2r*tmp2r + tmp2i*tmp2i;
2824 a = (tmp1r*tmp2r + tmp1i*tmp2i) / tmp;
2825 b = (tmp1i*tmp2r - tmp1r*tmp2i) / tmp;
2826 b2 = b*b;
2827 b = 2*a*b;
2828 a = a*a-b2;
2829
2830 a2 = Sqr(a);
2831 b2 = Sqr(b);
2832 dist2 = a2+b2;
2833
2834 if(dist2 < mindist2) mindist2 = dist2;
2835 tmp1r = a-1;
2836 if(dist2 > 10000.0 || tmp1r*tmp1r+b2 < 1/10000.0)
2837 {
2838 cf = fractal_exterior_color(TPat, col, a, b);
2839 break;
2840 }
2841 }
2842
2843 if(col == it_max)
2844 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2845
2846 return(cf);
2847 }
2848
2849
2850 /*****************************************************************************
2851 *
2852 * FUNCTION
2853 *
2854 * magnet2m_pattern
2855 *
2856 * INPUT
2857 *
2858 * EPoint -- The point in 3d space at which the pattern
2859 * is evaluated.
2860 *
2861 * OUTPUT
2862 *
2863 * RETURNS
2864 *
2865 * DBL value in the range 0.0 to 1.0
2866 *
2867 * AUTHOR
2868 *
2869 * Nieminen Juha
2870 *
2871 * DESCRIPTION
2872 *
2873 * -
2874 *
2875 * CHANGES
2876 *
2877 * -
2878 *
2879 ******************************************************************************/
2880
magnet2m_pattern(VECTOR EPoint,TPATTERN * TPat)2881 static DBL magnet2m_pattern (VECTOR EPoint, TPATTERN *TPat)
2882 {
2883 int it_max, col;
2884 DBL a, b, cf
2885 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2886 =0
2887 #endif
2888 , a2, b2, x, y, tmp, tmp1r, tmp1i, tmp2r, tmp2i,
2889 c1r, c2r, c1c2r, c1c2i, dist2, mindist2;
2890
2891 x = EPoint[X];
2892 y = EPoint[Y];
2893 a = a2 = 0;
2894 b = b2 = 0;
2895 mindist2 = 10000;
2896
2897 c1r = x-1; c2r = x-2;
2898 c1c2r = c1r*c2r-y*y;
2899 c1c2i = (c1r+c2r)*y;
2900
2901 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2902
2903 for (col = 0; col < it_max; col++)
2904 {
2905 tmp1r = a2*a-3*a*b2 + 3*(a*c1r-b*y) + c1c2r;
2906 tmp1i = 3*a2*b-b2*b + 3*(a*y+b*c1r) + c1c2i;
2907 tmp2r = 3*(a2-b2) + 3*(a*c2r-b*y) + c1c2r + 1;
2908 tmp2i = 6*a*b + 3*(a*y+b*c2r) + c1c2i;
2909 tmp = tmp2r*tmp2r + tmp2i*tmp2i;
2910 a = (tmp1r*tmp2r + tmp1i*tmp2i) / tmp;
2911 b = (tmp1i*tmp2r - tmp1r*tmp2i) / tmp;
2912 b2 = b*b;
2913 b = 2*a*b;
2914 a = a*a-b2;
2915
2916 a2 = Sqr(a);
2917 b2 = Sqr(b);
2918 dist2 = a2+b2;
2919
2920 if(dist2 < mindist2) mindist2 = dist2;
2921 tmp1r = a-1;
2922 if(dist2 > 10000.0 || tmp1r*tmp1r+b2 < 1/10000.0)
2923 {
2924 cf = fractal_exterior_color(TPat, col, a, b);
2925 break;
2926 }
2927 }
2928
2929 if(col == it_max)
2930 cf = fractal_interior_color(TPat, col, a, b, mindist2);
2931
2932 return(cf);
2933 }
2934
2935
2936 /*****************************************************************************
2937 *
2938 * FUNCTION
2939 *
2940 * magnet2j_pattern
2941 *
2942 * INPUT
2943 *
2944 * EPoint -- The point in 3d space at which the pattern
2945 * is evaluated.
2946 *
2947 * OUTPUT
2948 *
2949 * RETURNS
2950 *
2951 * DBL value in the range 0.0 to 1.0
2952 *
2953 * AUTHOR
2954 *
2955 * Nieminen Juha
2956 *
2957 * DESCRIPTION
2958 *
2959 * -
2960 *
2961 * CHANGES
2962 *
2963 * -
2964 *
2965 ******************************************************************************/
2966
magnet2j_pattern(VECTOR EPoint,TPATTERN * TPat)2967 static DBL magnet2j_pattern (VECTOR EPoint, TPATTERN *TPat)
2968 {
2969 int it_max, col;
2970 DBL a, b, cf
2971 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
2972 =0
2973 #endif
2974 , a2, b2, tmp, tmp1r, tmp1i, tmp2r, tmp2i, c1r,c2r,c1c2r,c1c2i,
2975 cr = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[U], ci = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Coord[V],
2976 dist2, mindist2;
2977
2978 a = EPoint[X]; a2 = Sqr(a);
2979 b = EPoint[Y]; b2 = Sqr(b);
2980 mindist2 = a2+b2;
2981
2982 c1r = cr-1, c2r = cr-2;
2983 c1c2r = c1r*c2r-ci*ci;
2984 c1c2i = (c1r+c2r)*ci;
2985
2986 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
2987
2988 for (col = 0; col < it_max; col++)
2989 {
2990 tmp1r = a2*a-3*a*b2 + 3*(a*c1r-b*ci) + c1c2r;
2991 tmp1i = 3*a2*b-b2*b + 3*(a*ci+b*c1r) + c1c2i;
2992 tmp2r = 3*(a2-b2) + 3*(a*c2r-b*ci) + c1c2r + 1;
2993 tmp2i = 6*a*b + 3*(a*ci+b*c2r) + c1c2i;
2994 tmp = tmp2r*tmp2r + tmp2i*tmp2i;
2995 a = (tmp1r*tmp2r + tmp1i*tmp2i) / tmp;
2996 b = (tmp1i*tmp2r - tmp1r*tmp2i) / tmp;
2997 b2 = b*b;
2998 b = 2*a*b;
2999 a = a*a-b2;
3000
3001 a2 = Sqr(a);
3002 b2 = Sqr(b);
3003 dist2 = a2+b2;
3004
3005 if(dist2 < mindist2) mindist2 = dist2;
3006 tmp1r = a-1;
3007 if(dist2 > 10000.0 || tmp1r*tmp1r+b2 < 1/10000.0)
3008 {
3009 cf = fractal_exterior_color(TPat, col, a, b);
3010 break;
3011 }
3012 }
3013
3014 if(col == it_max)
3015 cf = fractal_interior_color(TPat, col, a, b, mindist2);
3016
3017 return(cf);
3018 }
3019
3020
3021 /*****************************************************************************
3022 *
3023 * FUNCTION
3024 *
3025 * mandel_pattern
3026 *
3027 * INPUT
3028 *
3029 * EPoint -- The point in 3d space at which the pattern
3030 * is evaluated.
3031 *
3032 * OUTPUT
3033 *
3034 * RETURNS
3035 *
3036 * DBL value in the range 0.0 to 1.0
3037 *
3038 * AUTHOR
3039 *
3040 * submitted by user, name lost (sorry)
3041 *
3042 * DESCRIPTION
3043 *
3044 * The mandel pattern computes the standard Mandelbrot fractal pattern and
3045 * projects it onto the X-Y plane. It uses the X and Y coordinates to compute
3046 * the Mandelbrot set.
3047 *
3048 * CHANGES
3049 *
3050 * Oct 1994 : adapted from pigment by [CY]
3051 * May 2001 : updated with code from Warp [trf]
3052 *
3053 ******************************************************************************/
3054
mandel_pattern(VECTOR EPoint,TPATTERN * TPat)3055 static DBL mandel_pattern (VECTOR EPoint, TPATTERN *TPat)
3056 {
3057 int it_max, col;
3058 DBL a, b, cf
3059 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
3060 =0
3061 #endif
3062 , a2, b2, x, y, dist2, mindist2;
3063
3064 a = x = EPoint[X]; a2 = Sqr(a);
3065 b = y = EPoint[Y]; b2 = Sqr(b);
3066 mindist2 = a2+b2;
3067
3068 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
3069
3070 for (col = 0; col < it_max; col++)
3071 {
3072 b = 2.0 * a * b + y;
3073 a = a2 - b2 + x;
3074
3075 a2 = Sqr(a);
3076 b2 = Sqr(b);
3077 dist2 = a2+b2;
3078
3079 if(dist2 < mindist2) mindist2 = dist2;
3080 if(dist2 > 4.0)
3081 {
3082 cf = fractal_exterior_color(TPat, col, a, b);
3083 break;
3084 }
3085 }
3086
3087 if(col == it_max)
3088 cf = fractal_interior_color(TPat, col, a, b, mindist2);
3089
3090 return(cf);
3091 }
3092
3093
3094 /*****************************************************************************
3095 *
3096 * FUNCTION
3097 *
3098 * mandel3_pattern
3099 *
3100 * INPUT
3101 *
3102 * EPoint -- The point in 3d space at which the pattern
3103 * is evaluated.
3104 *
3105 * OUTPUT
3106 *
3107 * RETURNS
3108 *
3109 * DBL value in the range 0.0 to 1.0
3110 *
3111 * AUTHOR
3112 *
3113 * Nieminen Juha
3114 *
3115 * DESCRIPTION
3116 *
3117 * -
3118 *
3119 * CHANGES
3120 *
3121 * -
3122 *
3123 ******************************************************************************/
3124
mandel3_pattern(VECTOR EPoint,TPATTERN * TPat)3125 static DBL mandel3_pattern (VECTOR EPoint, TPATTERN *TPat)
3126 {
3127 int it_max, col;
3128 DBL a, b, cf
3129 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
3130 =0
3131 #endif
3132 , a2, b2, x, y, dist2, mindist2;
3133
3134 a = x = EPoint[X]; a2 = Sqr(a);
3135 b = y = EPoint[Y]; b2 = Sqr(b);
3136 mindist2 = a2+b2;
3137
3138 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
3139
3140 for (col = 0; col < it_max; col++)
3141 {
3142 b = 3.0*a2*b - b2*b + y;
3143 a = a2*a - 3.0*a*b2 + x;
3144
3145 a2 = Sqr(a);
3146 b2 = Sqr(b);
3147 dist2 = a2+b2;
3148
3149 if(dist2 < mindist2) mindist2 = dist2;
3150 if(dist2 > 4.0)
3151 {
3152 cf = fractal_exterior_color(TPat, col, a, b);
3153 break;
3154 }
3155 }
3156
3157 if(col == it_max)
3158 cf = fractal_interior_color(TPat, col, a, b, mindist2);
3159
3160 return(cf);
3161 }
3162
3163
3164 /*****************************************************************************
3165 *
3166 * FUNCTION
3167 *
3168 * mandel4_pattern
3169 *
3170 * INPUT
3171 *
3172 * EPoint -- The point in 3d space at which the pattern
3173 * is evaluated.
3174 *
3175 * OUTPUT
3176 *
3177 * RETURNS
3178 *
3179 * DBL value in the range 0.0 to 1.0
3180 *
3181 * AUTHOR
3182 *
3183 * Nieminen Juha
3184 *
3185 * DESCRIPTION
3186 *
3187 * -
3188 *
3189 * CHANGES
3190 *
3191 * -
3192 *
3193 ******************************************************************************/
3194
mandel4_pattern(VECTOR EPoint,TPATTERN * TPat)3195 static DBL mandel4_pattern (VECTOR EPoint, TPATTERN *TPat)
3196 {
3197 int it_max, col;
3198 DBL a, b, cf
3199 #ifdef AVOID_MIGHT_BE_USED_UNINITIALIZED_WARNINGS_PATCH
3200 =0
3201 #endif
3202 , a2, b2, x, y, dist2, mindist2;
3203
3204 a = x = EPoint[X]; a2 = Sqr(a);
3205 b = y = EPoint[Y]; b2 = Sqr(b);
3206 mindist2 = a2+b2;
3207
3208 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
3209
3210 for (col = 0; col < it_max; col++)
3211 {
3212 b = 4.0 * (a2*a*b - a*b2*b) + y;
3213 a = a2*a2 - 6.0*a2*b2 + b2*b2 + x;
3214
3215 a2 = Sqr(a);
3216 b2 = Sqr(b);
3217 dist2 = a2+b2;
3218
3219 if(dist2 < mindist2) mindist2 = dist2;
3220 if(dist2 > 4.0)
3221 {
3222 cf = fractal_exterior_color(TPat, col, a, b);
3223 break;
3224 }
3225 }
3226
3227 if(col == it_max)
3228 cf = fractal_interior_color(TPat, col, a, b, mindist2);
3229
3230 return(cf);
3231 }
3232
3233
3234 /*****************************************************************************
3235 *
3236 * FUNCTION
3237 *
3238 * mandelx_pattern
3239 *
3240 * INPUT
3241 *
3242 * EPoint -- The point in 3d space at which the pattern
3243 * is evaluated.
3244 *
3245 * OUTPUT
3246 *
3247 * RETURNS
3248 *
3249 * DBL value in the range 0.0 to 1.0
3250 *
3251 * AUTHOR
3252 *
3253 * Nieminen Juha
3254 *
3255 * DESCRIPTION
3256 *
3257 * -
3258 *
3259 * CHANGES
3260 *
3261 * -
3262 *
3263 ******************************************************************************/
3264
mandelx_pattern(VECTOR EPoint,TPATTERN * TPat)3265 static DBL mandelx_pattern (VECTOR EPoint, TPATTERN *TPat)
3266 {
3267 int it_max, col, exponent;
3268 DBL a, b, cf=0, x, y, dist2, mindist2;
3269 int* binomial_coeff;
3270
3271 if(BinomialCoefficientsInited == false)
3272 InitializeBinomialCoefficients();
3273
3274 a = x = EPoint[X];
3275 b = y = EPoint[Y];
3276 mindist2 = a*a+b*b;
3277
3278 it_max = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
3279 exponent = TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Exponent;
3280
3281 binomial_coeff = &BinomialCoefficients[(exponent+1)*exponent/2];
3282
3283 for (col = 0; col < it_max; col++)
3284 {
3285 // Calculate (a+bi)^exponent
3286 DBL new_a = pow(a, exponent);
3287 for(int k=2; k<=exponent; k+=2)
3288 {
3289 new_a += binomial_coeff[k]*pow(a, exponent-k)*pow(b, k);
3290 }
3291 DBL new_b = 0;
3292 for(int l=1; l<=exponent; l+=2)
3293 {
3294 new_b += binomial_coeff[l]*pow(a, exponent-l)*pow(b, l);
3295 }
3296
3297 a = new_a + x;
3298 b = new_b + y;
3299
3300 dist2 = a*a+b*b;
3301
3302 if(dist2 < mindist2) mindist2 = dist2;
3303 if(dist2 > 4.0)
3304 {
3305 cf = fractal_exterior_color(TPat, col, a, b);
3306 break;
3307 }
3308 }
3309
3310 if(col == it_max)
3311 cf = fractal_interior_color(TPat, col, a, b, mindist2);
3312
3313 return(cf);
3314 }
3315
3316
3317 /*****************************************************************************
3318 *
3319 * FUNCTION
3320 *
3321 * marble_pattern
3322 *
3323 * INPUT
3324 *
3325 * EPoint -- The point in 3d space at which the pattern
3326 * is evaluated.
3327 * TPat -- Texture pattern struct
3328 *
3329 * OUTPUT
3330 *
3331 * RETURNS
3332 *
3333 * DBL value in the range 0.0 to 1.0
3334 *
3335 * AUTHOR
3336 *
3337 * POV-Ray Team
3338 *
3339 * DESCRIPTION
3340 *
3341 * CHANGES
3342 *
3343 * Oct 1994 : adapted from pigment by [CY]
3344 *
3345 ******************************************************************************/
3346
marble_pattern(VECTOR EPoint,TPATTERN * TPat)3347 static DBL marble_pattern (VECTOR EPoint, TPATTERN *TPat)
3348 {
3349 register DBL turb_val;
3350 TURB *Turb;
3351
3352 if ((Turb=Search_For_Turb(TPat->Warps)) != NULL)
3353 {
3354 turb_val = Turb->Turbulence[X] * Turbulence(EPoint,Turb,TPat);
3355 }
3356 else
3357 {
3358 turb_val = 0.0;
3359 }
3360
3361 return(EPoint[X] + turb_val);
3362 }
3363
3364
3365 /*****************************************************************************
3366 *
3367 * FUNCTION
3368 *
3369 * object_pattern
3370 *
3371 * INPUT
3372 *
3373 * EPoint -- The point in 3d space at which the pattern
3374 * is evaluated.
3375 * TPat -- Texture pattern struct
3376 *
3377 * OUTPUT
3378 *
3379 * RETURNS
3380 *
3381 * DBL value in the range 0.0 to 1.0
3382 *
3383 * AUTHOR
3384 *
3385 * DESCRIPTION
3386 *
3387 * CHANGES
3388 *
3389 ******************************************************************************/
3390
object_pattern(VECTOR EPoint,TPATTERN * TPat)3391 static DBL object_pattern (VECTOR EPoint, TPATTERN *TPat)
3392 {
3393 if(TPat->Vals.Object != NULL)
3394 {
3395 if(Inside_Object(EPoint, TPat->Vals.Object))
3396 return 1.0;
3397 else
3398 return 0.0;
3399 }
3400
3401 return 0.0;
3402 }
3403
3404 /*****************************************************************************
3405 *
3406 * FUNCTION
3407 *
3408 * onion_pattern
3409 *
3410 * INPUT
3411 *
3412 * EPoint -- The point in 3d space at which the pattern
3413 * is evaluated.
3414 *
3415 * OUTPUT
3416 *
3417 * RETURNS
3418 *
3419 * DBL value in the range 0.0 to 1.0
3420 *
3421 * AUTHOR
3422 *
3423 * Scott Taylor
3424 *
3425 * DESCRIPTION
3426 *
3427 * CHANGES
3428 *
3429 * Jul 1991 : Creation.
3430 * Oct 1994 : adapted from pigment by [CY]
3431 *
3432 ******************************************************************************/
3433
onion_pattern(VECTOR EPoint)3434 static DBL onion_pattern (VECTOR EPoint)
3435 {
3436 /* The variable noise is not used as noise in this function */
3437
3438 register DBL noise;
3439
3440 /*
3441 This ramp goes 0-1,1-0,0-1,1-0...
3442
3443 noise = (fmod(sqrt(Sqr(x)+Sqr(y)+Sqr(z)),2.0)-1.0);
3444
3445 if (noise<0.0) {noise = 0.0-noise;}
3446 */
3447
3448 /* This ramp goes 0-1, 0-1, 0-1, 0-1 ... */
3449
3450 noise = (fmod(sqrt(Sqr(EPoint[X])+Sqr(EPoint[Y])+Sqr(EPoint[Z])), 1.0));
3451
3452 return(noise);
3453 }
3454
3455
3456 /*****************************************************************************
3457 *
3458 * FUNCTION
3459 *
3460 * INPUT
3461 *
3462 * OUTPUT
3463 *
3464 * RETURNS
3465 *
3466 * AUTHOR
3467 *
3468 * DESCRIPTION
3469 *
3470 * CHANGES
3471 *
3472 ******************************************************************************/
3473
pigment_pattern(VECTOR EPoint,TPATTERN * TPat,INTERSECTION * isect)3474 static DBL pigment_pattern (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *isect)
3475 {
3476 DBL value;
3477 COLOUR Col;
3478 int colour_found=false;
3479
3480 if (TPat->Vals.Pigment)
3481 colour_found = Compute_Pigment(Col, TPat->Vals.Pigment, EPoint, isect);
3482
3483 if(!colour_found)
3484 value = 0.0;
3485 else
3486 value = GREY_SCALE(Col);
3487
3488 return value ;
3489 }
3490
3491
3492 /*****************************************************************************
3493 *
3494 * FUNCTION
3495 *
3496 * planar_pattern
3497 *
3498 * INPUT
3499 *
3500 * EPoint -- The point in 3d space at which the pattern
3501 * is evaluated.
3502 *
3503 * OUTPUT
3504 *
3505 * RETURNS
3506 *
3507 * DBL value in the range 0.0 to 1.0
3508 *
3509 * AUTHOR
3510 *
3511 * -
3512 *
3513 * DESCRIPTION
3514 *
3515 * -
3516 *
3517 * CHANGES
3518 *
3519 * -
3520 *
3521 ******************************************************************************/
3522
planar_pattern(VECTOR EPoint)3523 static DBL planar_pattern (VECTOR EPoint)
3524 {
3525 register DBL value = fabs(EPoint[Y]);
3526
3527 CLIP_DENSITY(value);
3528
3529 return value;
3530 }
3531
3532
3533 /*****************************************************************************
3534 *
3535 * FUNCTION
3536 *
3537 * quilted_pattern
3538 *
3539 * INPUT
3540 *
3541 * OUTPUT
3542 *
3543 * RETURNS
3544 *
3545 * AUTHOR
3546 *
3547 * Dan Farmer & Chris Young
3548 *
3549 * DESCRIPTION
3550 *
3551 * CHANGES
3552 *
3553 ******************************************************************************/
3554
quilted_pattern(VECTOR EPoint,TPATTERN * TPat)3555 static DBL quilted_pattern (VECTOR EPoint, TPATTERN *TPat)
3556 {
3557 VECTOR value;
3558 DBL t;
3559
3560 value[X] = EPoint[X]-FLOOR(EPoint[X])-0.5;
3561 value[Y] = EPoint[Y]-FLOOR(EPoint[Y])-0.5;
3562 value[Z] = EPoint[Z]-FLOOR(EPoint[Z])-0.5;
3563
3564 t = sqrt(value[X]*value[X]+value[Y]*value[Y]+value[Z]*value[Z]);
3565
3566 t = quilt_cubic(t, TPat->Vals.Quilted LESS_MEMORY_PATCH_MO Control0, TPat->Vals.Quilted LESS_MEMORY_PATCH_MO Control1);
3567
3568 value[X] *= t;
3569 value[Y] *= t;
3570 value[Z] *= t;
3571
3572 return((fabs(value[X])+fabs(value[Y])+fabs(value[Z]))/3.0);
3573 }
3574
3575
3576 /*****************************************************************************
3577 *
3578 * FUNCTION
3579 *
3580 * radial_pattern
3581 *
3582 * INPUT
3583 *
3584 * EPoint -- The point in 3d space at which the pattern
3585 * is evaluated.
3586 *
3587 * OUTPUT
3588 *
3589 * RETURNS
3590 *
3591 * DBL value in the range 0.0 to 1.0
3592 *
3593 * AUTHOR
3594 *
3595 * Chris Young -- new in vers 2.0
3596 *
3597 * DESCRIPTION
3598 *
3599 * CHANGES
3600 *
3601 * Oct 1994 : adapted from pigment by [CY]
3602 *
3603 ******************************************************************************/
3604
radial_pattern(VECTOR EPoint)3605 static DBL radial_pattern (VECTOR EPoint)
3606 {
3607 register DBL value;
3608
3609 if ((fabs(EPoint[X])<0.001) && (fabs(EPoint[Z])<0.001))
3610 {
3611 value = 0.25;
3612 }
3613 else
3614 {
3615 value = 0.25 + (atan2(EPoint[X],EPoint[Z]) + M_PI) / TWO_M_PI;
3616 }
3617
3618 return(value);
3619 }
3620
3621
3622 /*****************************************************************************
3623 *
3624 * FUNCTION
3625 *
3626 * ripples_pattern
3627 *
3628 * INPUT
3629 *
3630 * EPoint -- The point in 3d space at which the pattern
3631 * is evaluated.
3632 * TPat -- Texture pattern struct
3633 *
3634 * OUTPUT
3635 *
3636 * RETURNS
3637 *
3638 * DBL value in the range 0.0 to 1.0
3639 *
3640 * AUTHOR
3641 *
3642 * POV-Ray Team
3643 *
3644 * DESCRIPTION : Note this pattern is only used for pigments and textures.
3645 * Normals have a specialized pattern for this.
3646 *
3647 * CHANGES
3648 *
3649 * Nov 1994 : adapted from normal by [CY]
3650 *
3651 ******************************************************************************/
3652
ripples_pattern(VECTOR EPoint,TPATTERN * TPat)3653 static DBL ripples_pattern (VECTOR EPoint, TPATTERN *TPat)
3654 {
3655 register unsigned int i;
3656 register DBL length, index;
3657 DBL scalar =0.0;
3658 VECTOR point;
3659
3660 for (i = 0 ; i < Number_Of_Waves ; i++)
3661 {
3662 VSub (point, EPoint, Wave_Sources[i]);
3663 VLength (length, point);
3664
3665 if (length == 0.0)
3666 length = 1.0;
3667
3668 index = length * TPat->Frequency + TPat->Phase;
3669
3670 scalar += cycloidal(index);
3671 }
3672
3673 scalar = 0.5*(1.0+(scalar / (DBL)Number_Of_Waves));
3674
3675 return(scalar);
3676 }
3677
3678
3679 /*****************************************************************************
3680 *
3681 * FUNCTION
3682 *
3683 * slope_pattern
3684 *
3685 * INPUT
3686 *
3687 * EPoint -- The point in 3d space at which the pattern
3688 * is evaluated.
3689 * TPat -- Texture pattern struct
3690 * Intersection - intersection struct
3691 *
3692 * OUTPUT
3693 *
3694 * RETURNS
3695 *
3696 * DBL value in the range 0.0 to 1.0, 0.0 if normal is NULL
3697 *
3698 * AUTHOR
3699 *
3700 * -hdf-
3701 *
3702 * DESCRIPTION :
3703 *
3704 * calculates the surface slope from surface normal vector
3705 *
3706 * CHANGES
3707 *
3708 * Apr 1998 : written by H.-D. Fink
3709 * May 1998 : modified by M.C. Andrews - now combines slope and 'gradient'.
3710 *
3711 ******************************************************************************/
3712
slope_pattern(VECTOR EPoint,TPATTERN * TPat,INTERSECTION * Isection)3713 static DBL slope_pattern (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *Isection)
3714 {
3715 DBL value, value1, value2;
3716
3717 if (Isection == NULL) return 0.0; /* just in case ... */
3718
3719 if (TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Base > 0)
3720 /* short case 1: slope vector in x, y or z direction */
3721 value1 = Isection->PNormal[TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Base - 1];
3722 else if (TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Base < 0)
3723 /* short case 2: slope vector in negative x, y or z direction */
3724 value1 = -Isection->PNormal[-TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Base - 1];
3725 else
3726 /* projection slope onto normal vector */
3727 VDot(value1, Isection->PNormal, TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Vector);
3728
3729 /* Clamp to 1.0. */
3730 /* should never be necessary since both vectors are normalized */
3731 if (value1 > 1.0) value1 = 1.0;
3732 else if (value1 < -1.0) value1 = -1.0;
3733
3734 value1 = asin(value1) / M_PI * 2;
3735 value1 = (value1 + 1.0) * 0.5; /* normalize to [0..1] interval */
3736
3737 /* If set, use offset and scalings for slope and altitude. */
3738 if (0.0 != TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Mod[V])
3739 {
3740 value1 = (value1 - TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Mod[U]) / TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Mod[V];
3741 }
3742
3743 if (!TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Len)
3744 {
3745 /* Clamp to 1.0. */
3746 if ( value1 == 1.0 )
3747 {
3748 value1= value1- EPSILON;
3749 }
3750 else
3751 {
3752 value1 = (value1 < 0.0) ? 1.0 + fmod(value1, 1.0) : fmod(value1, 1.0);
3753 }
3754 return value1; /* no altitude defined */
3755 }
3756
3757 /* Calculate projection of Epoint along altitude vector */
3758 if (TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Base > 0)
3759 /* short case 1: altitude vector in x, y or z direction */
3760 value2 = EPoint[TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Base - 1];
3761 else if (TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Base < 0)
3762 /* short case 2: altitude vector in negative x, y or z direction */
3763 value2 = -EPoint[-TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Base - 1];
3764 else
3765 /* projection of Epoint along altitude vector */
3766 VDot(value2, EPoint, TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Vector);
3767
3768 if (0.0 != TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Mod[V])
3769 {
3770 value2 = (value2 - TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Mod[U]) / TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Mod[V];
3771 }
3772
3773 value = TPat->Vals.Slope LESS_MEMORY_PATCH_MO Slope_Len * value1 + TPat->Vals.Slope LESS_MEMORY_PATCH_MO Altit_Len * value2;
3774
3775 /* Clamp to 1.0. */
3776 if ( value - 1.0 < EPSILON && value >= 1.0 )
3777 {
3778 /* 1.0 is a very common value to get *exactly*. We don't want to wrap
3779 it to the bottom end of the map. */
3780 value = value - EPSILON;
3781 }
3782 else
3783 {
3784 value = (value < 0.0) ? 1.0 + fmod(value, 1.0) : fmod(value, 1.0);
3785 }
3786 return value;
3787
3788 }
3789
3790
3791 /*****************************************************************************
3792 *
3793 * FUNCTION
3794 *
3795 * spiral1_pattern
3796 *
3797 * INPUT
3798 *
3799 * EPoint -- The point in 3d space at which the pattern
3800 * is evaluated.
3801 * TPat -- Texture pattern struct
3802 *
3803 * OUTPUT
3804 *
3805 * RETURNS
3806 *
3807 * DBL value in the range 0.0 to 1.0
3808 *
3809 * AUTHOR
3810 *
3811 * Dieter Bayer
3812 *
3813 * DESCRIPTION
3814 * Spiral whirles around z-axis.
3815 * The number of "arms" is defined in the TPat.
3816 *
3817 * CHANGES
3818 *
3819 * Aug 1994 : Creation.
3820 * Oct 1994 : adapted from pigment by [CY]
3821 *
3822 ******************************************************************************/
3823
spiral1_pattern(VECTOR EPoint,TPATTERN * TPat)3824 static DBL spiral1_pattern (VECTOR EPoint, TPATTERN *TPat)
3825 {
3826 DBL rad, phi, turb_val;
3827 DBL x = EPoint[X];
3828 DBL y = EPoint[Y];
3829 DBL z = EPoint[Z];
3830 TURB *Turb;
3831
3832 if ((Turb=Search_For_Turb(TPat->Warps)) != NULL)
3833 {
3834 turb_val = Turb->Turbulence[X] * Turbulence(EPoint,Turb,TPat);
3835 }
3836 else
3837 {
3838 turb_val = 0.0;
3839 }
3840
3841 /* Get distance from z-axis. */
3842
3843 rad = sqrt(x * x + y * y);
3844
3845 /* Get angle in x,y-plane (0...2 PI). */
3846
3847 if (rad == 0.0)
3848 {
3849 phi = 0.0;
3850 }
3851 else
3852 {
3853 if (x < 0.0)
3854 {
3855 phi = 3.0 * M_PI_2 - asin(y / rad);
3856 }
3857 else
3858 {
3859 phi = M_PI_2 + asin(y / rad);
3860 }
3861 }
3862
3863 return(z + rad + (DBL)TPat->Vals.Arms * phi / TWO_M_PI + turb_val);
3864 }
3865
3866
3867 /*****************************************************************************
3868 *
3869 * FUNCTION
3870 *
3871 * spiral2_pattern
3872 *
3873 * INPUT
3874 *
3875 * EPoint -- The point in 3d space at which the pattern
3876 * is evaluated.
3877 * TPat -- Texture pattern struct
3878 *
3879 * OUTPUT
3880 *
3881 * RETURNS
3882 *
3883 * DBL value in the range 0.0 to 1.0
3884 *
3885 * AUTHOR
3886 *
3887 * Dieter Bayer
3888 *
3889 * DESCRIPTION
3890 * Spiral whirles around z-axis.
3891 * The number of "arms" is defined in the TPat.
3892 *
3893 * CHANGES
3894 *
3895 * Aug 1994 : Creation.
3896 * Oct 1994 : adapted from pigment by [CY]
3897 *
3898 ******************************************************************************/
3899
spiral2_pattern(VECTOR EPoint,TPATTERN * TPat)3900 static DBL spiral2_pattern (VECTOR EPoint, TPATTERN *TPat)
3901 {
3902 DBL rad, phi, turb_val;
3903 DBL x = EPoint[X];
3904 DBL y = EPoint[Y];
3905 DBL z = EPoint[Z];
3906 TURB *Turb;
3907
3908 if ((Turb=Search_For_Turb(TPat->Warps)) != NULL)
3909 {
3910 turb_val = Turb->Turbulence[X] * Turbulence(EPoint,Turb,TPat);
3911 }
3912 else
3913 {
3914 turb_val = 0.0;
3915 }
3916
3917 /* Get distance from z-axis. */
3918
3919 rad = sqrt(x * x + y * y);
3920
3921 /* Get angle in x,y-plane (0...2 PI) */
3922
3923 if (rad == 0.0)
3924 {
3925 phi = 0.0;
3926 }
3927 else
3928 {
3929 if (x < 0.0)
3930 {
3931 phi = 3.0 * M_PI_2 - asin(y / rad);
3932 }
3933 else
3934 {
3935 phi = M_PI_2 + asin(y / rad);
3936 }
3937 }
3938
3939 turb_val = Triangle_Wave(z + rad + (DBL)TPat->Vals.Arms * phi / TWO_M_PI +
3940 turb_val);
3941
3942 return(Triangle_Wave(rad) + turb_val);
3943 }
3944
3945
3946 /*****************************************************************************
3947 *
3948 * FUNCTION
3949 *
3950 * spherical_pattern
3951 *
3952 * INPUT
3953 *
3954 * EPoint -- The point in 3d space at which the pattern
3955 * is evaluated.
3956 *
3957 * OUTPUT
3958 *
3959 * RETURNS
3960 *
3961 * DBL value in the range 0.0 to 1.0
3962 *
3963 * AUTHOR
3964 *
3965 * -
3966 *
3967 * DESCRIPTION
3968 *
3969 * -
3970 *
3971 * CHANGES
3972 *
3973 * -
3974 *
3975 ******************************************************************************/
3976
spherical_pattern(VECTOR EPoint)3977 static DBL spherical_pattern (VECTOR EPoint)
3978 {
3979 register DBL value;
3980
3981 VLength(value, EPoint);
3982 CLIP_DENSITY(value);
3983
3984 return(value);
3985 }
3986
3987
3988 /*****************************************************************************
3989 *
3990 * FUNCTION
3991 *
3992 * waves_pattern
3993 *
3994 * INPUT
3995 *
3996 * EPoint -- The point in 3d space at which the pattern
3997 * is evaluated.
3998 * TPat -- Texture pattern struct
3999 *
4000 * OUTPUT
4001 *
4002 * RETURNS
4003 *
4004 * DBL value in the range 0.0 to 1.0
4005 *
4006 * AUTHOR
4007 *
4008 * POV-Ray Team
4009 *
4010 * DESCRIPTION : Note this pattern is only used for pigments and textures.
4011 * Normals have a specialized pattern for this.
4012 *
4013 * CHANGES
4014 *
4015 * Nov 1994 : adapted from normal by [CY]
4016 *
4017 ******************************************************************************/
4018
waves_pattern(VECTOR EPoint,TPATTERN * TPat)4019 static DBL waves_pattern (VECTOR EPoint, TPATTERN *TPat)
4020 {
4021 register unsigned int i;
4022 register DBL length, index;
4023 DBL scalar = 0.0;
4024 VECTOR point;
4025
4026 for (i = 0 ; i < Number_Of_Waves ; i++)
4027 {
4028 VSub (point, EPoint, Wave_Sources[i]);
4029 VLength (length, point);
4030
4031 if (length == 0.0)
4032 {
4033 length = 1.0;
4034 }
4035
4036 index = length * TPat->Frequency * frequency[i] + TPat->Phase;
4037
4038 scalar += cycloidal(index)/frequency[i];
4039 }
4040
4041 scalar = 0.2*(2.5+(scalar / (DBL)Number_Of_Waves));
4042
4043 return(scalar);
4044 }
4045
4046
4047 /*****************************************************************************
4048 *
4049 * FUNCTION
4050 *
4051 * wood_pattern
4052 *
4053 * INPUT
4054 *
4055 * EPoint -- The point in 3d space at which the pattern
4056 * is evaluated.
4057 * TPat -- Texture pattern struct
4058 *
4059 * OUTPUT
4060 *
4061 * RETURNS
4062 *
4063 * DBL value in the range 0.0 to 1.0
4064 *
4065 * AUTHOR
4066 *
4067 * POV-Ray Team
4068 *
4069 * DESCRIPTION
4070 *
4071 * CHANGES
4072 *
4073 * Oct 1994 : adapted from pigment by [CY]
4074 *
4075 ******************************************************************************/
4076
wood_pattern(VECTOR EPoint,TPATTERN * TPat)4077 static DBL wood_pattern (VECTOR EPoint, TPATTERN *TPat)
4078 {
4079 register DBL length;
4080 VECTOR WoodTurbulence;
4081 VECTOR point;
4082 DBL x=EPoint[X];
4083 DBL y=EPoint[Y];
4084 TURB *Turb;
4085
4086 if ((Turb=Search_For_Turb(TPat->Warps)) != NULL)
4087 {
4088 DTurbulence (WoodTurbulence, EPoint,Turb);
4089 point[X] = cycloidal((x + WoodTurbulence[X]) * Turb->Turbulence[X]);
4090 point[Y] = cycloidal((y + WoodTurbulence[Y]) * Turb->Turbulence[Y]);
4091 }
4092 else
4093 {
4094 point[X] = 0.0;
4095 point[Y] = 0.0;
4096 }
4097 point[Z] = 0.0;
4098
4099 point[X] += x;
4100 point[Y] += y;
4101
4102 /* point[Z] += z; Deleted per David Buck -- BP 7/91 */
4103
4104 VLength (length, point);
4105
4106 return(length);
4107 }
4108
4109
4110 /*****************************************************************************
4111 *
4112 * FUNCTION
4113 *
4114 * wrinkles_pattern
4115 *
4116 * INPUT
4117 *
4118 * EPoint -- The point in 3d space at which the pattern
4119 * is evaluated.
4120 *
4121 * OUTPUT
4122 *
4123 * RETURNS
4124 *
4125 * DBL value in the range 0.0 to 1.0
4126 *
4127 * AUTHOR
4128 *
4129 * POV-Ray Team
4130 *
4131 * DESCRIPTION : Note this pattern is only used for pigments and textures.
4132 * Normals have a specialized pattern for this.
4133 *
4134 * CHANGES
4135 *
4136 * Nov 1994 : adapted from normal by [CY]
4137 *
4138 ******************************************************************************/
4139
wrinkles_pattern(VECTOR EPoint,TPATTERN * TPat)4140 static DBL wrinkles_pattern (VECTOR EPoint, TPATTERN *TPat)
4141 {
4142 register int i;
4143 DBL lambda = 2.0;
4144 DBL omega = 0.5;
4145 DBL value;
4146 VECTOR temp;
4147 DBL noise;
4148 int noise_generator = 0;
4149
4150 if (TPat != NULL)
4151 noise_generator = (TPat->Flags & NOISE_FLAGS) >> 4;
4152 if (noise_generator == 0)
4153 noise_generator = opts.Noise_Generator;
4154
4155 if(noise_generator>1)
4156 {
4157 noise = Noise(EPoint, TPat)*2.0-0.5;
4158 value = min(max(noise,0.0),1.0);
4159 }
4160 else
4161 {
4162 value = Noise(EPoint, TPat);
4163 }
4164
4165 for (i = 1; i < 10; i++)
4166 {
4167 VScale(temp,EPoint,lambda);
4168
4169 if(noise_generator>1)
4170 {
4171 noise = Noise(temp, TPat)*2.0-0.5;
4172 value += omega * min(max(noise,0.0),1.0);
4173 }
4174 else
4175 {
4176 value += omega * Noise(temp, TPat);
4177 }
4178
4179 lambda *= 2.0;
4180
4181 omega *= 0.5;
4182 }
4183
4184 return(value/2.0);
4185 }
4186
4187
4188 /*****************************************************************************
4189 *
4190 * FUNCTION
4191 *
4192 * IntPickInCube(tvx,tvy,tvz, p1)
4193 * a version of PickInCube that takes integers for input
4194 *
4195 * INPUT
4196 *
4197 * ?
4198 *
4199 * OUTPUT
4200 *
4201 * RETURNS
4202 *
4203 * long integer hash function used, to speed up cacheing.
4204 *
4205 * AUTHOR
4206 *
4207 * original PickInCube by Jim McElhiney
4208 * this integer one modified by Nathan Kopp
4209 *
4210 * DESCRIPTION
4211 *
4212 * A subroutine to go with crackle.
4213 *
4214 * Pick a random point in the same unit-sized cube as tv, in a
4215 * predictable way, so that when called again with another point in
4216 * the same unit cube, p1 is picked to be the same.
4217 *
4218 * CHANGES
4219 *
4220 ******************************************************************************/
4221
IntPickInCube(int tvx,int tvy,int tvz,VECTOR p1)4222 static int IntPickInCube(int tvx, int tvy, int tvz, VECTOR p1)
4223 {
4224 int seed, temp;
4225
4226 seed = Hash3d(tvx&0xFFF,tvy&0xFFF,tvz&0xFFF);
4227
4228 temp = POV_GET_OLD_RAND(); /* save current seed */
4229 POV_SRAND(seed);
4230
4231 p1[X] = tvx + FRAND();
4232 p1[Y] = tvy + FRAND();
4233 p1[Z] = tvz + FRAND();
4234
4235 POV_SRAND(temp); /* restore */
4236
4237 return((int)seed);
4238 }
4239
4240
4241 /*****************************************************************************
4242 *
4243 * FUNCTION
4244 *
4245 * PickInCube(tv, p1)
4246 *
4247 * INPUT
4248 *
4249 * ?
4250 *
4251 * OUTPUT
4252 *
4253 * RETURNS
4254 *
4255 * long integer hash function used, to speed up cacheing.
4256 *
4257 * AUTHOR
4258 *
4259 * Jim McElhiney
4260 *
4261 * DESCRIPTION
4262 *
4263 * A subroutine to go with crackle.
4264 *
4265 * Pick a random point in the same unit-sized cube as tv, in a
4266 * predictable way, so that when called again with another point in
4267 * the same unit cube, p1 is picked to be the same.
4268 *
4269 * CHANGES
4270 *
4271 ******************************************************************************/
4272
PickInCube(VECTOR tv,VECTOR p1)4273 int PickInCube(VECTOR tv, VECTOR p1)
4274 {
4275 int seed, temp;
4276 VECTOR flo;
4277
4278 /*
4279 * This uses floor() not FLOOR, so it will not be a mirror
4280 * image about zero in the range -1.0 to 1.0. The viewer
4281 * won't see an artefact around the origin.
4282 */
4283
4284 flo[X] = floor(tv[X] - EPSILON);
4285 flo[Y] = floor(tv[Y] - EPSILON);
4286 flo[Z] = floor(tv[Z] - EPSILON);
4287
4288 seed = Hash3d((int)flo[X], (int)flo[Y], (int)flo[Z]);
4289
4290 temp = POV_GET_OLD_RAND(); /* save current seed */
4291
4292 POV_SRAND(seed);
4293
4294 p1[X] = flo[X] + FRAND();
4295 p1[Y] = flo[Y] + FRAND();
4296 p1[Z] = flo[Z] + FRAND();
4297
4298 POV_SRAND(temp); /* restore */
4299
4300 return((int)seed);
4301 }
4302
4303
4304 /*****************************************************************************
4305 *
4306 * FUNCTION
4307 *
4308 * INPUT
4309 *
4310 * OUTPUT
4311 *
4312 * RETURNS
4313 *
4314 * AUTHOR
4315 *
4316 * DESCRIPTION
4317 *
4318 * CHANGES
4319 *
4320 ******************************************************************************/
4321
4322 const DBL INV_SQRT_3_4 = 1.154700538;
quilt_cubic(DBL t,DBL p1,DBL p2)4323 DBL quilt_cubic(DBL t, DBL p1, DBL p2)
4324 {
4325 DBL it=(1-t);
4326 DBL itsqrd=it*it;
4327 /* DBL itcubed=it*itsqrd; */
4328 DBL tsqrd=t*t;
4329 DBL tcubed=t*tsqrd;
4330 DBL val;
4331
4332 /* Originally coded as...
4333
4334 val= (DBL)(itcubed*n1+(tcubed)*n2+3*t*(itsqrd)*p1+3*(tsqrd)*(it)*p2);
4335
4336 re-written by CEY to optimise because n1=0 n2=1 always.
4337
4338 */
4339
4340 val = (tcubed + 3.0*t*itsqrd*p1 + 3.0*tsqrd*it*p2) * INV_SQRT_3_4;
4341
4342 return(val);
4343 }
4344
4345
4346 /*****************************************************************************
4347 *
4348 * FUNCTION
4349 *
4350 * INPUT
4351 *
4352 * OUTPUT
4353 *
4354 * RETURNS
4355 *
4356 * AUTHOR
4357 *
4358 * DESCRIPTION
4359 *
4360 * CHANGES
4361 *
4362 ******************************************************************************/
4363
fractal_exterior_color(TPATTERN * TPat,int iters,DBL a,DBL b)4364 static DBL fractal_exterior_color(TPATTERN *TPat, int iters, DBL a, DBL b)
4365 {
4366 switch(TPat->Vals.Fractal LESS_MEMORY_PATCH_MO exterior_type)
4367 {
4368 case 0:
4369 return (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4370 case 1:
4371 return (DBL)iters / (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO Iterations;
4372 case 2:
4373 return a * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4374 case 3:
4375 return b * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4376 case 4:
4377 return a*a * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4378 case 5:
4379 return b*b * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4380 case 6:
4381 return sqrt(a*a+b*b) * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO efactor;
4382 }
4383 return 0;
4384 }
4385
4386
4387 /*****************************************************************************
4388 *
4389 * FUNCTION
4390 *
4391 * INPUT
4392 *
4393 * OUTPUT
4394 *
4395 * RETURNS
4396 *
4397 * AUTHOR
4398 *
4399 * DESCRIPTION
4400 *
4401 * CHANGES
4402 *
4403 ******************************************************************************/
4404
fractal_interior_color(TPATTERN * TPat,int iters,DBL a,DBL b,DBL mindist2)4405 static DBL fractal_interior_color(TPATTERN *TPat, int iters, DBL a, DBL b, DBL mindist2)
4406 {
4407 switch(TPat->Vals.Fractal LESS_MEMORY_PATCH_MO interior_type)
4408 {
4409 case 0:
4410 return (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4411 case 1:
4412 return sqrt(mindist2) * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4413 case 2:
4414 return a * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4415 case 3:
4416 return b * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4417 case 4:
4418 return a*a * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4419 case 5:
4420 return b*b * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4421 case 6:
4422 return a*a+b*b * (DBL)TPat->Vals.Fractal LESS_MEMORY_PATCH_MO ifactor;
4423 }
4424 return 0;
4425 }
4426
4427
4428 /*****************************************************************************
4429 *
4430 * FUNCTION
4431 *
4432 * Create_Density_File
4433 *
4434 * INPUT
4435 *
4436 * OUTPUT
4437 *
4438 * RETURNS
4439 *
4440 * AUTHOR
4441 *
4442 * Dieter Bayer
4443 *
4444 * DESCRIPTION
4445 *
4446 * Create a density file structure.
4447 *
4448 * CHANGES
4449 *
4450 * Dec 1996 : Creation.
4451 *
4452 ******************************************************************************/
4453
Create_Density_File()4454 DENSITY_FILE *Create_Density_File()
4455 {
4456 DENSITY_FILE *New;
4457
4458 New = (DENSITY_FILE *)POV_MALLOC(sizeof(DENSITY_FILE), "density file");
4459
4460 New->Interpolation = NO_INTERPOLATION;
4461
4462 New->Data = (DENSITY_FILE_DATA *)POV_MALLOC(sizeof(DENSITY_FILE_DATA), "density file data");
4463
4464 New->Data->References = 1;
4465
4466 New->Data->Name = NULL;
4467
4468 New->Data->Sx =
4469 New->Data->Sy =
4470 New->Data->Sz = 0;
4471
4472 New->Data->Type = 0;
4473
4474 New->Data->Density32 = NULL;
4475 New->Data->Density16 = NULL;
4476 New->Data->Density8 = NULL;
4477
4478 return (New);
4479 }
4480
4481
4482 /*****************************************************************************
4483 *
4484 * FUNCTION
4485 *
4486 * Copy_Density_File
4487 *
4488 * INPUT
4489 *
4490 * OUTPUT
4491 *
4492 * RETURNS
4493 *
4494 * AUTHOR
4495 *
4496 * Dieter Bayer
4497 *
4498 * DESCRIPTION
4499 *
4500 * Copy a density file structure.
4501 *
4502 * CHANGES
4503 *
4504 * Dec 1996 : Creation.
4505 *
4506 ******************************************************************************/
4507
Copy_Density_File(DENSITY_FILE * Old)4508 DENSITY_FILE *Copy_Density_File(DENSITY_FILE *Old)
4509 {
4510 DENSITY_FILE *New;
4511
4512 if (Old != NULL)
4513 {
4514 New = (DENSITY_FILE *)POV_MALLOC(sizeof(DENSITY_FILE), "density file");
4515
4516 *New = *Old;
4517
4518 New->Data->References++;
4519 }
4520 else
4521 {
4522 New=NULL;
4523 }
4524
4525 return(New);
4526 }
4527
4528
4529 /*****************************************************************************
4530 *
4531 * FUNCTION
4532 *
4533 * Destroy_Density_File
4534 *
4535 * INPUT
4536 *
4537 * OUTPUT
4538 *
4539 * RETURNS
4540 *
4541 * AUTHOR
4542 *
4543 * Dieter Bayer
4544 *
4545 * DESCRIPTION
4546 *
4547 * Destroy a density file structure.
4548 *
4549 * CHANGES
4550 *
4551 * Dec 1996 : Creation.
4552 *
4553 ******************************************************************************/
4554
Destroy_Density_File(DENSITY_FILE * Density_File)4555 void Destroy_Density_File(DENSITY_FILE *Density_File)
4556 {
4557 if(Density_File != NULL)
4558 {
4559 if((--(Density_File->Data->References)) == 0)
4560 {
4561 POV_FREE(Density_File->Data->Name);
4562
4563 if(Density_File->Data->Type == 4)
4564 {
4565 POV_FREE(Density_File->Data->Density32);
4566 }
4567 else if(Density_File->Data->Type == 2)
4568 {
4569 POV_FREE(Density_File->Data->Density16);
4570 }
4571 else if(Density_File->Data->Type == 1)
4572 {
4573 POV_FREE(Density_File->Data->Density8);
4574 }
4575
4576 POV_FREE(Density_File->Data);
4577 }
4578
4579 POV_FREE(Density_File);
4580 }
4581 }
4582
Read_Density_File(DENSITY_FILE * df)4583 void Read_Density_File(DENSITY_FILE *df)
4584 {
4585 int x, y, z, sx, sy, sz;
4586 IStream *file;
4587 unsigned long len;
4588
4589 if (df == NULL)
4590 return;
4591
4592 /* Allocate and read density file. */
4593
4594 if((df != NULL) && (df->Data->Name != NULL))
4595 {
4596 if((file = Locate_File(df->Data->Name, POV_File_Data_DF3, NULL, true)) == NULL)
4597 Error("Cannot read media density file.");
4598
4599 sx = df->Data->Sx = readushort(file);
4600 sy = df->Data->Sy = readushort(file);
4601 sz = df->Data->Sz = readushort(file);
4602
4603 file->seekg(0, IOBase::seek_end);
4604 len = file->tellg() - 6;
4605 file->seekg(6);
4606
4607 // figure out the data size
4608 #ifdef AVOID_TYPE_CONVERSION_WARNINGS_PATCH
4609 if(((unsigned long)sx * sy * sz * 4) == len)
4610 #else
4611 if((sx * sy * sz * 4) == len)
4612 #endif
4613 {
4614 df->Data->Type = 4;
4615
4616 unsigned int *map = (unsigned int *)POV_MALLOC(sx * sy * sz * sizeof(unsigned int), "media density file data 32 bit");
4617
4618 for (z = 0; z < sz; z++)
4619 {
4620 for (y = 0; y < sy; y++)
4621 {
4622 for (x = 0; x < sx; x++)
4623 map[z * sy * sx + y * sx + x] = readuint(file);
4624 }
4625 }
4626
4627 df->Data->Density32 = map;
4628 }
4629 #ifdef AVOID_TYPE_CONVERSION_WARNINGS_PATCH
4630 else if(((unsigned long)sx * sy * sz * 2) == len)
4631 #else
4632 else if((sx * sy * sz * 2) == len)
4633 #endif
4634 {
4635 df->Data->Type = 2;
4636
4637 unsigned short *map = (unsigned short *)POV_MALLOC(sx * sy * sz * sizeof(unsigned short), "media density file data 16 bit");
4638
4639 for (z = 0; z < sz; z++)
4640 {
4641 for (y = 0; y < sy; y++)
4642 {
4643 for (x = 0; x < sx; x++)
4644 map[z * sy * sx + y * sx + x] = readushort(file);
4645 }
4646 }
4647
4648 df->Data->Density16 = map;
4649 }
4650 #ifdef AVOID_TYPE_CONVERSION_WARNINGS_PATCH
4651 else if(((unsigned long)sx * sy * sz) == len)
4652 #else
4653 else if((sx * sy * sz) == len)
4654 #endif
4655 {
4656 df->Data->Type = 1;
4657
4658 unsigned char *map = (unsigned char *)POV_MALLOC(sx * sy * sz * sizeof(unsigned char), "media density file data 8 bit");
4659
4660 for (z = 0; z < sz; z++)
4661 {
4662 for (y = 0; y < sy; y++)
4663 file->read((char *)(&(map[z * sy * sx + y * sx])), sizeof(unsigned char) * sx);
4664 }
4665
4666 df->Data->Density8 = map;
4667 }
4668 else
4669 Error("Invalid density file size");
4670
4671 if (file != NULL)
4672 {
4673 delete file;
4674 }
4675 }
4676 }
4677
readushort(IStream * infile)4678 static unsigned short readushort(IStream *infile)
4679 {
4680 short i0 = 0, i1 = 0;
4681
4682 #ifdef AVOID_ALWAYS_FALSE_COMPARISON_WARNINGS_PATCH
4683 i0 = infile->Read_Byte ();
4684 i1 = infile->Read_Byte ();
4685 #else
4686 if ((i0 = infile->Read_Byte ()) == EOF || (i1 = infile->Read_Byte ()) == EOF)
4687 {
4688 Error("Error reading density_file");
4689 }
4690 #endif
4691
4692 return (((unsigned short)i0 << 8) | (unsigned short)i1);
4693 }
4694
readuint(IStream * infile)4695 static unsigned int readuint(IStream *infile)
4696 {
4697 int i0 = 0, i1 = 0, i2 = 0, i3 = 0;
4698
4699 #ifdef AVOID_ALWAYS_FALSE_COMPARISON_WARNINGS_PATCH
4700 i0 = infile->Read_Byte ();
4701 i1 = infile->Read_Byte ();
4702 i2 = infile->Read_Byte ();
4703 i3 = infile->Read_Byte ();
4704 #else
4705 if ((i0 = infile->Read_Byte ()) == EOF || (i1 = infile->Read_Byte ()) == EOF ||
4706 (i2 = infile->Read_Byte ()) == EOF || (i3 = infile->Read_Byte ()) == EOF)
4707 {
4708 Error("Error reading density_file");
4709 }
4710 #endif
4711
4712 return (((unsigned int)i0 << 24) | ((unsigned int)i1 << 16) | ((unsigned int)i2 << 8) | (unsigned int)i3);
4713 }
4714
4715 // This should be called once, either at povray start or the first time
4716 // a fractal pattern is created:
InitializeBinomialCoefficients()4717 static void InitializeBinomialCoefficients()
4718 {
4719 int* ptr = BinomialCoefficients;
4720 *ptr = 1; ++ptr;
4721
4722 for(unsigned n=1; n<=FRACTAL_MAX_EXPONENT; ++n)
4723 {
4724 *ptr = 1; ++ptr;
4725 for(unsigned k=1; k<n; ++k)
4726 {
4727 *ptr = *(ptr-(n+1)) + *(ptr-n); ++ptr;
4728 }
4729 *ptr = 1; ++ptr;
4730 }
4731 ptr = BinomialCoefficients+1;
4732 for(unsigned m=1; m<=FRACTAL_MAX_EXPONENT; ++m)
4733 {
4734 ++ptr;
4735 for(unsigned k=1; k<m; ++k)
4736 {
4737 if((k&2)!=0) *ptr = -(*ptr);
4738 ++ptr;
4739 }
4740 if((m&2)!=0) *ptr = -(*ptr);
4741 ++ptr;
4742 }
4743
4744 BinomialCoefficientsInited = true;
4745 }
4746
4747
4748 //patches from here on
4749 #ifdef LISTED_PATTERN_PATCH
ListedPat(VECTOR EPoint,TPATTERN * TPat)4750 static DBL ListedPat(VECTOR EPoint, TPATTERN *TPat)
4751 {
4752 return TPat->Vals.ListedVal;
4753 }
4754 #endif
4755
4756 /*****************************************************************************
4757 *
4758 * FUNCTION
4759 *
4760 * angle_of_incidence
4761 *
4762 * INPUT
4763 *
4764 * EPoint -- The point in 3d space at which the pattern
4765 * is evaluated.
4766 * TPat -- Texture pattern struct
4767 * Intersection - intersection struct
4768 *
4769 * OUTPUT
4770 *
4771 * RETURNS
4772 *
4773 * DBL value in the range 0.0 to 1.0
4774 *
4775 * AUTHOR
4776 *
4777 * ML
4778 *
4779 * DESCRIPTION :
4780 *
4781 * calculates the angle of incidence between surface normal vector and a point
4782 *
4783 * CHANGES
4784 * Feb 2001 : written by ML
4785 *
4786 ******************************************************************************/
4787
4788 #ifdef ANGLE_OF_INCIDENCE_PATCH
4789
angle_of_incidence(VECTOR EPoint,TPATTERN * TPat,INTERSECTION * Intersection)4790 static DBL angle_of_incidence (VECTOR EPoint, TPATTERN *TPat, INTERSECTION *Intersection)
4791 {
4792 DBL value;
4793 VECTOR v_origin; /* vector EPoint to origin */
4794
4795 if (Intersection == NULL) return 0.0; /* just in case ... */
4796
4797 if (TPat->Vals.Aoi LESS_MEMORY_PATCH_MO pt_fixed)
4798 {
4799 /* angle / fixed point */
4800 VSub(v_origin,TPat->Vals.Aoi LESS_MEMORY_PATCH_MO AOI_origin,EPoint);
4801 if ((v_origin[X]!=0) || (v_origin[Y]!=0) || (v_origin[Z]!=0))
4802 VNormalizeEq(v_origin);
4803 }
4804 else
4805 {
4806 /* angle / ray direction */
4807 v_origin[X]=-(Intersection->RDir)[X];
4808 v_origin[Y]=-(Intersection->RDir)[Y];
4809 v_origin[Z]=-(Intersection->RDir)[Z];
4810 }
4811
4812 /* Cos(Angle_of_Incidence) */
4813 VDot(value, Intersection->PNormal, v_origin);
4814
4815 /* Angle from cos and normalize to 0..1 interval */
4816 value=acos(value)/M_PI;
4817
4818 return value;
4819
4820 }
4821 #endif
4822
4823 /*****************************************************************************
4824 *
4825 * FUNCTION
4826 *
4827 * projection_pat
4828 *
4829 * INPUT
4830 *
4831 * EPoint -- The point in 3d space at which the pattern
4832 * is evaluated.
4833 * TPat -- Texture pattern struct
4834 *
4835 * OUTPUT
4836 *
4837 * RETURNS
4838 *
4839 * DBL value 0.0 or 1.0
4840 *
4841 * AUTHOR
4842 *
4843 * ML
4844 *
4845 * DESCRIPTION :
4846 *
4847 * trace a ray from EPoint in direction dirv,
4848 * return 1 if hit object 0 else
4849 *
4850 * CHANGES
4851 * Feb 2001 : written by ML
4852 *
4853 ******************************************************************************/
4854
4855 #ifdef PROJECTION_PATCH
4856
projection_pat(VECTOR EPoint,TPATTERN * TPat,INTERSECTION * Inter)4857 static DBL projection_pat (VECTOR EPoint, TPATTERN *TPat,INTERSECTION *Inter)
4858 {
4859 INTERSECTION Intersect;
4860 RAY Ray;
4861 VECTOR vsave,Jitter_Offset;
4862 DBL x,y,z,r,t;
4863 DBL value;
4864 int i;
4865
4866
4867 Initialize_Ray_Containers( &Ray );
4868 Assign_Vector( Ray.Initial, EPoint );
4869
4870 switch (TPat->Vals.projection LESS_MEMORY_PATCH_MO type)
4871 {
4872 case 0:
4873 VSub(Ray.Direction,TPat->Vals.projection LESS_MEMORY_PATCH_MO dirv,EPoint);
4874 /* test if not null vector */
4875 if ( (Ray.Direction[X]==0) && (Ray.Direction[Y]==0) && (Ray.Direction[Z]==0) )
4876 return 0;
4877 VNormalizeEq(Ray.Direction);
4878 break;
4879 case 1:
4880 Assign_Vector( Ray.Direction, TPat->Vals.projection LESS_MEMORY_PATCH_MO dirv);
4881 break;
4882 case 2:
4883 if (!Inter) return 0.0;
4884 Assign_Vector( Ray.Direction, Inter->PNormal);
4885 break;
4886 }
4887
4888 if (TPat->Vals.projection LESS_MEMORY_PATCH_MO samples)
4889 {
4890 /* blur */
4891 value=0;
4892
4893 /* save our vector */
4894 Assign_Vector( vsave, Ray.Direction);
4895
4896 /* sum for all samples */
4897 for (i=0;i<TPat->Vals.projection LESS_MEMORY_PATCH_MO samples;i++)
4898 {
4899 /* Pick a random point on the surface of a unit sphere.
4900 Uses a method called the "trig method" explained in the
4901 comp.graphics.algorithms FAQ. */
4902 z = (FRAND() * 2) - 1;
4903 t = FRAND() * M_PI * 2;
4904 r = sqrt(1 - z*z);
4905 x = r * cos(t);
4906 y = r * sin(t);
4907 Make_Vector(Jitter_Offset, x, y, z);
4908
4909 /* Add the jittering to the vector */
4910 Assign_Vector(Ray.Direction, vsave);
4911 VAddScaledEq(Ray.Direction,TPat->Vals.projection LESS_MEMORY_PATCH_MO blurr, Jitter_Offset);
4912
4913 /* trace */
4914 if ( Intersection( &Intersect, TPat->Vals.projection LESS_MEMORY_PATCH_MO obj , &Ray ) )
4915 value++;
4916 }
4917 /* average */
4918 value/=TPat->Vals.projection LESS_MEMORY_PATCH_MO samples;
4919 }
4920 else
4921 {
4922 /* no blur */
4923 if ( Intersection( &Intersect, TPat->Vals.projection LESS_MEMORY_PATCH_MO obj , &Ray ) )
4924 value=1;
4925 else
4926 value=0;
4927 }
4928
4929 return value;
4930 }
4931 #endif
4932 END_POV_NAMESPACE
4933