Rendering 3D fractals without distance estimators

[Syntopia]

  Nice render.  For whatever reason, I'm having a few problems with the coloring function- you got it to do what I wanted it to do.  I'll mess around with it and take into consideration your post above this one. 

  I had a few problems with atan as well, trying to get atan2 type outputs.

  I think I might have that figured out- had to use some math trick, but am not satisfied using the math trick as the GLSL documentation site says that the atan function should output atan2 type outputs (I'll go back and retry atan(y,x) and atan (y/x) to see if they are different as the documentation seems to imply). 

  I'm assuming that GPU rendering with trigonometric functions is pretty fast , even on my old beast. 

In GLSL, the two-argument atan(y,x) is an atan2 type function. For x>0 it is equal to the one-argument atan(y/x).

I've attached my whole code, in case you want to see it.

Code:

// Output generated from file: Unnamed
// Created: ti 13. nov 21:32:37 2012

#info Mandelbulb without Distance Estimator

#define providesInside

#define providesColor

#include "Brute-Raytracer.frag"
#group Mandelbulb
//#define IterationsBetweenRedraws 5      // I really need to clean this code up... LATER though...

// Number of fractal iterations.
uniform int Iterations;  slider[0,4,100]

// Number of color iterations.
uniform int ColorIterations;  slider[0,3,100]
uniform int ColorIterationsMin;  slider[0,1,100]

// Mandelbulb exponent (8 is standard)
//uniform float Power; slider[-10,8,10]

// Bailout radius
uniform float Bailout; slider[0,12,30]

//uniform vec3 RotVector; slider[(0,0,0),(1,1,1),(1,1,1)]

//uniform float RotAngle; slider[0.00,0,180]

uniform bool Julia; checkbox[false]
uniform vec3 JuliaC; slider[(-2,-2,-2),(0,0,0),(2,2,2)]

vec3 cycle(vec3 c, float s) {
	return vec3(0.5)+0.5*vec3(cos(s*Cycles+c.x),cos(s*Cycles+c.y),cos(s*Cycles+c.z));
}

vec3 color(vec3 p) {
	orbitTrap = vec4(10000.0);
	inside(p);

	orbitTrap.w = sqrt(orbitTrap.w);
	vec3 orbitColor;
	if (CycleColors) {
		orbitColor = cycle(X.xyz,orbitTrap.x)*X.w*orbitTrap.x +
		cycle(Y.xyz,orbitTrap.y)*Y.w*orbitTrap.y +
		cycle(Z.xyz,orbitTrap.z)*Z.w*orbitTrap.z +
		cycle(R.xyz,orbitTrap.w)*R.w*orbitTrap.w;
	} else {
		orbitColor = X.xyz*X.w*orbitTrap.x +
		Y.xyz*Y.w*orbitTrap.y +
		Z.xyz*Z.w*orbitTrap.z +
		R.xyz*R.w*orbitTrap.w;
	}
	
	//vec3 color = mix(BaseColor, 3.0*orbitColor,  OrbitStrength);
	return orbitColor;
	
}
void powN2(inout vec3 z, float zr0) {
	float sr23=sqrt(2./3.);
	float sr13=sqrt(1./3.);
	float nx=z.x*sr23-z.z*sr13;
	float sz=z.x*sr13 + z.z*sr23;
	float sx=nx;
	float sr12=sqrt(.5);
	nx=sx*sr12-z.y*sr12;             
	float sy=sx*sr12+z.y*sr12;
	sx=nx*nx;
	sy=sy*sy;
	float ny=sy;
	sz=sz*sz;
	float r2=sx+sy+sz;
	if (r2!=0.) {                                       
		nx=(sx+r2)*(9.*sx-sy-sz)/(9.*sx+sy+sz)-.5;
		ny=(sy+r2)*(9.*sy-sx-sz)/(9.*sy+sx+sz)-.5;
		sz=(sz+r2)*(9.*sz-sx-sy)/(9.*sz+sx+sy)-.5;
		
	}
	sx=nx;
	sy=ny;
	nx=sx*sr12+sy*sr12;
	sy=-sx*sr12+sy*sr12; 
	sx=nx;
	nx=sx*sr23+sz*sr13;
	sz=-sx*sr13+sz*sr23;				//some things can be cleaned up
	sx=nx;
	float sx2=sx*sx;
	float sy2=sy*sy;				// will be switching code around later       
	float sz2=sz*sz;
	nx=sx2-sy2-sz2;
	float r3=2.*abs(sx)/sqrt(sy2+sz2);
	float nz=r3*(sy2-sz2);
 	ny=r3*2.*sy*sz;
	z= vec3(nx,ny,nz);
}


bool inside(vec3 pos) {
	vec3 z=pos;
	float r;
	int i=0;
	r=length(z);
	while(r<Bailout && (i<Iterations)) {
		powN2(z,r);
		if (i>ColorIterationsMin && i<ColorIterations) {
			orbitTrap=min(orbitTrap,vec4(z.x*z.x,z.y*z.y,z.z*z.z,dot(z,z)));

		}
		if (pos.x>0.) {                          // can I do this???
			z+=vec3(pos.x*.5,0.,0.);
		} else {
			z+=vec3(pos.x,0.,0.);
		}	
		r=length(z);
		
		i++;
	}
	return (r<Bailout);
}




#preset default
FOV = 0.62536
Eye = -2.06349,0.0288363,0.0261386
Target = 6.7978,-0.0399766,-0.0109628
Up = -0.0063931,-0.829538,0.0116395
EquiRectangular = false
Gamma = 2.5
Exposure = 1.34694
Brightness = 1
Contrast = 0.9901
Saturation = 1
SpecularExp = 5.455
SpotLightDir = 0.87654,0.78126
CamLight = 1,1,1,1.13044
CamLightMin = 0
BackgroundColor = 0.666667,0.666667,0.498039
ToneMapping = 3
Near = 0.12372
Far = 2.45904
NormalScale = 0.00024
AOScale = 0.00631
Glow = 0.34167
AOStrength = 0.86047
Samples = 100
Stratify = true
DebugInside = false
SampleNeighbors = true
Specular = 10
SpotLight = 1,0.678431,0.494118,0.78431
Fog = 1.84
ShowDepth = false
DebugNormals = false
BaseColor = 1,1,1
OrbitStrength = 1
X = 1,1,1,0.86408
Y = 0.345098,0.666667,0,1
Z = 1,0.666667,0,1
R = 0.0784314,1,0.941176,0.54902
GradientBackground = 0.76085
CycleColors = false
Cycles = 6.47914
Iterations = 10
ColorIterations = 3
Bailout = 12
Julia = false
JuliaC = 0,0,0
ColorIterationsMin = 1
#endpreset

[richardrosenman]

Hi Syntopia;

Thank you for that explanation regarding the color mapping to orbit. I'm excited to try it out.

I decided to whip up a comparison for your debugging needs - in case it helps:



These are both rendered with 1000 subframes (you should consider upping the limit of 2000 for the next release) with identical settings and  they have not been downsampled. However, I think the forum may be doing some of that (downsampling) so here's the original image:

http://www.richardrosenman.com/storage/brute-raytracer_comparison.jpg

Hope this helps!

-Rich

P.S. Darn. Of course your code above doesn't work with my new 'old' fragment files...

[M Benesi]

In GLSL, the two-argument atan(y,x) is an atan2 type function. For x>0 it is equal to the one-argument atan(y/x).

  Thanks.  Will check that out- must have been using atan(y/x) instead of atan(y,x).

I've attached my whole code, in case you want to see it.

  Cool.  Definitely will try it out.  Glad you saw to grab the orbit prior to adding in the x pixel component. 

[M Benesi]

  Ok.  Your program is awesome.  Seems likes it's more than an order of magnitude faster than CPU rendering.  Wow...   Also, if you want this split from this thread, I can delete it, a mod could move it, or whatever.  

  Anyways, here is some slightly improved code for <everyone> to play with, after a few bits of info.  I think the trig version is a bit slower than the algebraic, but I'm not sure.  If it is, I've noticed it takes half the time to use complex numbers in ChaosPro for integer "n", although I don't know how much faster that would be on the GPU.   ALSO-  Don't know if I should assign more variables in the following code to break it up into more streams- whoever is familiar with this stuff could say something (or maybe I could simply do a quick re-write and test it, ehh??)

  For the "Spokes" option set it to an integer value (otherwise you get breaks).  4 if you're going for a fake holy grail type fractal, 3 is normal.      Side note: couldn't get intBitsToFloat to work, so had to use a float instead of an integer for Spokes.

  There are also manPower and magPower sliders- these are for the 2 inter-related formulas, each with its own "n" (z^n).

  Note that the Mandy portion of the formula can be replaced with ANY of the various Mandy type formulas (like the 3d Burning Ship, the "flat" Mandy, etc.  Also, you can do lots of tricks with these formulas, like adding in the old portions from the beginning of a section to the end of the section, etc.  

  Lastly (well, not the ultimate lastly...) you can switch the Mandy portion to the beginning of the formula if you so desire.  In fact, you can separate the Mag portion with its rotations into a whole separate function call, then call various Mandy functions...

Code:

// Copied and modified from Syntopia's code   
// Created: ~ 10:30 pm 11-13-2012

#info Mandelbulb without Distance Estimator

#define providesInside

#define providesColor

#include "Brute-Raytracer.frag"
#group Mandelbulb
//#define IterationsBetweenRedraws 5      // I really need to clean this code up... LATER though...

// Number of fractal iterations.
uniform int Iterations;  slider[0,4,100]

// Number of color iterations.
uniform int ColorIterations;  slider[0,3,100]
uniform int ColorIterationsMin;  slider[0,1,100]

// Mandelbulb exponent (2 is standard)
uniform float manPower; slider[-10,2,10]

//mag exponent  (2 is standard)
uniform float magPower; slider[-10,2,10]

//  Spokes number 
uniform float Spokes; slider[-4,3,20]

// Bailout radius
uniform float Bailout; slider[0,12,30]

//uniform vec3 RotVector; slider[(0,0,0),(1,1,1),(1,1,1)]

//uniform float RotAngle; slider[0.00,0,180]

uniform bool Julia; checkbox[false]
uniform vec3 JuliaC; slider[(-2,-2,-2),(0,0,0),(2,2,2)]

vec3 cycle(vec3 c, float s) {
	return vec3(0.5)+0.5*vec3(cos(s*Cycles+c.x),cos(s*Cycles+c.y),cos(s*Cycles+c.z));
}

vec3 color(vec3 p) {
	orbitTrap = vec4(1.0);
	inside(p);

	orbitTrap.w = sqrt(orbitTrap.w);
	vec3 orbitColor;
	if (CycleColors) {
		orbitColor = cycle(X.xyz,orbitTrap.x)*X.w*orbitTrap.x +
		cycle(Y.xyz,orbitTrap.y)*Y.w*orbitTrap.y +
		cycle(Z.xyz,orbitTrap.z)*Z.w*orbitTrap.z +
		cycle(R.xyz,orbitTrap.w)*R.w*orbitTrap.w;
	} else {
		orbitColor = X.xyz*X.w*orbitTrap.x +
		Y.xyz*Y.w*orbitTrap.y +
		Z.xyz*Z.w*orbitTrap.z +
		R.xyz*R.w*orbitTrap.w;
	}
	
	//vec3 color = mix(BaseColor, 3.0*orbitColor,  OrbitStrength);
	return orbitColor;
	
}
void powN2(inout vec3 z, in float d) {
	float sr23=sqrt(2./3.);
	float sr13=sqrt(1./3.);
	float nx=z.x*sr23-z.z*sr13;
	float sz=z.x*sr13 + z.z*sr23;   // sz rotated
	float sx=nx;
	float sr12=sqrt(.5);
	nx=sx*sr12-z.y*sr12;               //nx
	float sy=sx*sr12+z.y*sr12;  //sy rotated
	sx=nx;

	float rxyz=pow((sx*sx+sy*sy+sz*sz),magPower*.5);

	float r1=sqrt(sy*sy+sz*sz);	
	float victor=atan(r1,abs(sx*d))*magPower;		// Is it faster to use
	nx=(sx*sx+rxyz)*cos(victor)-.5;				// multiple variables 

	r1=sqrt(sx*sx+sz*sz);						// for these to split
	victor=atan(r1,abs(sy*d))*magPower;		// processes???
	float ny=(sy*sy+rxyz)*cos(victor)-.5;

	r1=sqrt(sx*sx+sy*sy);
	victor=atan(r1,abs(sz*d))*magPower;
	sz=(sz*sz+rxyz)*cos(victor)-.5;

	sx=nx;
	sy=ny;

	nx=sx*sr12+sy*sr12;
	sy=-sx*sr12+sy*sr12; 
	sx=nx;
	nx=sx*sr23+sz*sr13;
	sz=-sx*sr13+sz*sr23;				//some things can be cleaned up
	sx=nx;

	rxyz=pow((sx*sx+sy*sy+sz*sz),manPower/2.);
	r1=sqrt(sy*sy+sz*sz);
	victor=atan(r1,sx)*manPower;
	float phi=atan(sz,sy)*manPower;
	sx=rxyz*cos(victor);
	r1=rxyz*sin(victor);
	sz=r1*cos(phi);
	sy=r1*sin(phi);
	z= vec3(sx,sy,sz);
	
}


bool inside(vec3 pos) {
	
		float ryz=sqrt(pos.y*pos.y+pos.z*pos.z);
		float theta=(Spokes)*atan(pos.z,pos.y)/3.;
		vec3 z=vec3(pos.x,ryz*cos(theta),ryz*sin(theta));
	 
	float r;
	int i=0;
	r=length(z);
	float d=magPower*1.5;						// sometimes I make d= magPower +1
	while(r<Bailout && (i<Iterations)) {			//  for a sharper fractal.. we can add a switch
		powN2(z,d);
		if (i>ColorIterationsMin && i<ColorIterations) {
			orbitTrap=min(orbitTrap,vec4(z.x*z.x,z.y*z.y,z.z*z.z,dot(z,z)));

		}
		if (pos.x>0.) {                          // can I do this???
			z+=vec3(pos.x*.5,0.,0.);
		} else {
			z+=vec3(pos.x,0.,0.);
		}	
		r=length(z);
		
		i++;
	}
	return (r<Bailout);
}




#preset default
FOV = 0.62536
Eye = 2.02493,-0.000258894,0.00365772
Target = -0.498369,0.00627817,-0.0886997
Up = -0.0024414,-0.999997,0
EquiRectangular = false
Gamma = 2.5
ToneMapping = 3
Exposure = 1.34694
Brightness = 1
Contrast = 0.9901
Saturation = 1
Near = 0.48984
Far = 2.71704
NormalScale = 0.00024
AOScale = 0.00291
Glow = 0.34167
AOStrength = 0.86047
Samples = 254
Stratify = true
DebugInside = false
SampleNeighbors = true
Specular = 1
SpecularExp = 5.455
SpotLight = 0.027451,0.196078,0.658824,0.89706
SpotLightDir = 0.80246,0.78126
CamLight = 1,1,1,1.44928
CamLightMin = 0.54217
Fog = 1.84
ShowDepth = false
DebugNormals = false
BaseColor = 0.388235,0.388235,0.388235
OrbitStrength = 0.78
X = 0.329412,0.321569,0.364706,0.08738
Y = 0.262745,0.333333,0.282353,0.02912
Z = 0.054902,0.0431373,0.439216,0.1068
R = 0.317647,0.670588,0.690196,0.64706
BackgroundColor = 0.0196078,0.0117647,0.192157
GradientBackground = 1.73915
CycleColors = false
Cycles = 23.1868
Iterations = 6
ColorIterations = 4
Bailout = 12
Julia = false
JuliaC = 0,0,0
ColorIterationsMin = 2
#endpreset

  Here are 2 images.  Click to enbiggen:

[Syntopia]

Hi Syntopia;

Thank you for that explanation regarding the color mapping to orbit. I'm excited to try it out.

I decided to whip up a comparison for your debugging needs - in case it helps:

Hi Richard,

I've changed from the hardware derivate version back to the original texture look up version. This removed a lot of jaggies. So sad - the hardware derivatives very much more elegant :-)

I've attached the modified DepthBufferShader.frag - please try it, and see if it improves the rendering. This also removes some of the artifacts when doing tile rendering - now it is only the SSAO, that makes visible tile artifacts.

<Quoted Image Removed>
These are both rendered with 1000 subframes (you should consider upping the limit of 2000 for the next release) with identical settings and  they have not been downsampled. However, I think the forum may be doing some of that (downsampling) so here's the original image:

Notice, that you could also increase the number of samples per subframe, from the default value of 1000.
Also notice, that if you go to continious mode, and put a zero in the Max (subframe) field in the main window, it will run indefinately.

[Syntopia]

 Ok.  Your program is awesome.  Seems likes it's more than an order of magnitude faster than CPU rendering.  Wow...   Also, if you want this split from this thread, I can delete it, a mod could move it, or whatever.  

Thanks, Benesi!

ALSO-  Don't know if I should assign more variables in the following code to break it up into more streams- whoever is familiar with this stuff could say something (or maybe I could simply do a quick re-write and test it, ehh??)

Not sure what you are asking, but use as few variables as posible, and try to use Vec3 wherever possible (the latter because of readability).
There is no need for the GPU to reorder instructions, or to try to make independent code paths. Each core of the GPU will execute all of your program (and they will run in lockstep in groups of 16 'threads' - executing the same command at the same time).

 For the "Spokes" option set it to an integer value (otherwise you get breaks).  4 if you're going for a fake holy grail type fractal, 3 is normal.      Side note: couldn't get intBitsToFloat to work, so had to use a float instead of an integer for Spokes.

To convert from integer to float, you just cast the numbers:

Code:

uniform int Spokes;  slider[0,9,100]

// Use Spokes as float
float f = float(Spokes)*2.0;

[richardrosenman]

Hi Syntopia;

Yes, the new fix seems to work correctly again -  nice!

I did some more rendering tests and I've noticed a couple of other things (should it be important).



I've been trying to get high resolution brute renderings despite the limitations. I've noticed rendering with tiles yields different lighting results than the same without tile rendering. Is this a typical limitation of tile rendering? In the example above, the first image was rendered at 4k, then downsampled to 1k. The second image was rendered at 1k without tile rendering. These do not have specular or occlusion enabled.

I also noticed that there is still tiling artefacts at the edges, even without spec or AO. It's not really visible in this example because it's been downsampled but I can post the original if required. You said there should be no more artefacts if not using AO? Here's the 4k original if you want to see it closely: Click here for 4k boring render.

Have a look and tell me what you think. I know you've mentioned the brute-renderer has many limitations when rendering in high-res but I thought I'd point these out for simply discussion's sake.

In full-screen rendering your new script works perfectly.

Link: http://www.richardrosenman.com/storage/tile_lighting_comparison.jpg

-Rich

[Syntopia]

I've been trying to get high resolution brute renderings despite the limitations. I've noticed rendering with tiles yields different lighting results than the same without tile rendering. Is this a typical limitation of tile rendering? In the example above, the first image was rendered at 4k, then downsampled to 1k. The second image was rendered at 1k without tile rendering. These do not have specular or occlusion enabled.

No, this is not typical - it must be an error. Looking at the pictures, I think the normals are somehow scaled scaled wrongly in the 'into screen' direction. The left image looks flattened.

I also noticed that there is still tiling artefacts at the edges, even without spec or AO. It's not really visible in this example because it's been downsampled but I can post the original if required. You said there should be no more artefacts if not using AO? Here's the 4k original if you want to see it closely: Click here for 4k boring render.

Yep, I said it, but I was wrong: The diffuse and specular part both depend on the normal, and the normal is also screen space calculated, so near the edges these will be wrong. You need to render only with ambient light, I guess :-(

I can fix the first part, but the second problem is a bit tougher. Still, it would be nice to have large scale renderings.

[richardrosenman]

I can fix the first part, but the second problem is a bit tougher. Still, it would be nice to have large scale renderings.

Yes, it would be but you're pretty damn close I'd say. The current tiling artefacts are minimal at best and easy Photoshop cleanup. The AO issue - well, no way around that I suppose if it uses screen-space AO unless a different approach to the AO is used. But I intend to do a high resolution final when this is all said and done so I will know how much Photoshop cleanup the AO requires.

I'd say the current big thing is the different lighting in tile rendering as you pointed out. But you sound confident you can fix that.

So really, it's pretty awesome as it is.

-Rich

[cKleinhuis]

if i shall split the arbitrary/double floating point topic, tell me from which posting to which

awesome work dudes, gotta play with it as well, and hack and slay in some formulas, lets get a nice hybrid renderer inside, cant wait to try it out, this will be featured in the last news of the year keep it coming

good nite, sleep well oO

[M Benesi]

Not sure what you are asking, but use as few variables as posible, and try to use Vec3 wherever possible (the latter because of readability).
There is no need for the GPU to reorder instructions, or to try to make independent code paths. Each core of the GPU will execute all of your program (and they will run in lockstep in groups of 16 'threads' - executing the same command at the same time).

  Italicized portion is the answer I was looking for.  Probably missed it in your blog.. or read it and forgot. 

To convert from integer to float, you just cast the numbers:

Code:

uniform int Spokes;  slider[0,9,100]

// Use Spokes as float
float f = float(Spokes)*2.0;

  Thanks- will do.  I'm forcing myself to take a break from math and code tonight, so it will probably be tomorrow (or later if something comes up).

[visual.bermarte]

@Benesi: great job!    - no distance estimation?

[M Benesi]

Thanks.  No DE type, as of yet.  Mikael's new brute force non-DE patch (from this thread) allowed this implementation.  I should look into DE for this type, simply because it would be fun to learn and useful for rendering. 

  I'm going to split off the one section of code so that people can see which section can be mixed with other codes- it might be useful as an intermediary function in KIFS- and I've some ideas to make it symmetric over all axes that I need to try out. 

[M Benesi]

To convert from integer to float, you just cast the numbers:

Code:

uniform int Spokes;  slider[0,9,100]

// Use Spokes as float
float f = float(Spokes)*2.0;

  For whatever reason (perhaps keyboard actuator failure on my part), the above code didn't work.  No problem however- not necessary for the fractal to work.

  Below this is code for a version with several rotations built in.  Note that using rotation 1 or 3 creates "breaks" in the fractal if you set spokes to anything other than the default value of 3.  Rotation 2 is AWESOME and works with various spokes values.  You should try it out... turns the whole fractal into a machine of fractal gears.  For Rotation 1 or 3, if you set the vector portion to 0,1,0 or any singular value like that, you get cool results as well.

 If I can figure out how to do animations with this awesome software, I'll throw together a quick demo- but my GPU is soooooo old....

  Here is the code- remember, Spokes setting other than 3 has discontinuities if you play with rotations 1 or 3.  Rotation 2 is awesome.

UPDATE Nov. 17th 2012:  If you set the magPower to 0, implemented a special mode that allows pure Mandy type.  Also added in a "seed" parameter- note that the code below I haven't tested (just adding these ideas off the top of my head- they work in ChaosPro).  Will update with modified TESTED code later.

Code:

// Output generated from file: Unnamed
// Created: ti 13. nov 21:32:37 2012

#info Mandelbulb without Distance Estimator

#define providesInside

#define providesColor
#include "MathUtils.frag"
#include "Brute-Raytracer.frag"
#group Fractal
//#define IterationsBetweenRedraws 5      // I really need to clean this code up... LATER though...

// Number of fractal iterations.
uniform int Iterations;  slider[0,4,100]

// Number of color iterations.
uniform int ColorIterations;  slider[0,3,100]
uniform int ColorIterationsMin;  slider[0,1,100]

// Mandelbulb exponent (2 is standard)
uniform float manPower; slider[-10,2,10]

//mag exponent  (2 is standard)
uniform float magPower; slider[-10,2,10]

//mag seed  (-.5 is standard)
uniform float magSeed; slider[-6,-.5,6]
//mag Softness
uniform float magSoft; slider[-10,1,100]

// Bailout radius
uniform float Bailout; slider[0,12,30]

uniform bool Julia; checkbox[false]
uniform vec3 JuliaC; slider[(-5,-5,-5),(0,0,0),(5,5,5)]

#group Rotation
//  Spokes number 
uniform float Spokes; slider[-4,3,20]

uniform vec3 RotVector1; slider[(0,0,0),(1,1,1),(1,1,1)]
uniform float RotAngle1; slider[-180.00,0,180]

mat3 rot1= rotationMatrix3(normalize(RotVector1), RotAngle1);

uniform vec3 RotVector2; slider[(0,0,0),(1,1,1),(1,1,1)]
uniform float RotAngle2; slider[-180.00,0,180]

mat3 rot2= rotationMatrix3(normalize(RotVector2), RotAngle2);

uniform vec3 RotVector3; slider[(0,0,0),(1,1,1),(1,1,1)]
uniform float RotAngle3; slider[-180.00,0,180]

mat3 rot3= rotationMatrix3(normalize(RotVector3), RotAngle3);


vec3 cycle(vec3 c, float s) {
	return vec3(0.5)+0.5*vec3(cos(s*Cycles+c.x),cos(s*Cycles+c.y),cos(s*Cycles+c.z));
}

vec3 color(vec3 p) {
	orbitTrap = vec4(1.0);
	inside(p);

	orbitTrap.w = sqrt(orbitTrap.w);
	vec3 orbitColor;
	if (CycleColors) {
		orbitColor = cycle(X.xyz,orbitTrap.x)*X.w*orbitTrap.x +
		cycle(Y.xyz,orbitTrap.y)*Y.w*orbitTrap.y +
		cycle(Z.xyz,orbitTrap.z)*Z.w*orbitTrap.z +
		cycle(R.xyz,orbitTrap.w)*R.w*orbitTrap.w;
	} else {
		orbitColor = X.xyz*X.w*orbitTrap.x +
		Y.xyz*Y.w*orbitTrap.y +
		Z.xyz*Z.w*orbitTrap.z +
		R.xyz*R.w*orbitTrap.w;
	}
	
	//vec3 color = mix(BaseColor, 3.0*orbitColor,  OrbitStrength);
	return orbitColor;
	
}
void magpart(inout vec3 z, in float d) {
					// from here until the next section break is the section that can
				 // be blended with other rotation based fractals- and maybe KIFS!

	float sr23=sqrt(2./3.);				
	float sr13=sqrt(1./3.);
	float nx=z.x*sr23-z.z*sr13;
	float sz=z.x*sr13 + z.z*sr23;   // sz rotated
	float sx=nx;
	float sr12=sqrt(.5);
	nx=sx*sr12-z.y*sr12;               //nx
	float sy=sx*sr12+z.y*sr12;  //sy rotated
	sx=nx;
	float ny=sy;
	float r1=0.;
	float rxyz=0.;
	float victor=0.;
	if ((magPower==2.)) {
		sx=nx*nx;
		sy=sy*sy;
	
		sz=sz*sz;
		r1=sx+sy+sz;
		if (r1!=0.) {                                       
			nx=(sx+r1)*(9.*sx-sy-sz)/(9.*sx+sy+sz)+magSeed;	
			ny=(sy+r1)*(9.*sy-sx-sz)/(9.*sy+sx+sz)+magSeed;	
			sz=(sz+r1)*(9.*sz-sx-sy)/(9.*sz+sx+sy)+magSeed;	
		}
	} else if ((magPower==0.)) {
		nx=(abs(sx)+magSeed)*2.;
		ny=(abs(sy)+magSeed)*2.;
		sz=(abs(sz)+magSeed)*2.;
	} else {
		rxyz=pow((sx*sx+sy*sy+sz*sz),(magPower)*.5);

		r1=sqrt(sy*sy+sz*sz);	
		victor=atan(r1,abs(sx*d))*magPower;		
		nx=(sx*sx+rxyz)*cos(victor)+magSeed;				//  The -.5 value is a seed value

		r1=sqrt(sx*sx+sz*sz);						//  you can use different ones
		victor=atan(r1,abs(sy*d))*magPower;		//  but ~ -.5 is the best, and keeping it the same for
		ny=(sy*sy+rxyz)*cos(victor)+magSeed;	//  all 3 variables works the best- but you can vary it between them

		r1=sqrt(sx*sx+sy*sy);
		victor=atan(r1,abs(sz*d))*magPower;
		sz=(sz*sz+rxyz)*cos(victor)+magSeed;	
	}
	vec3 z2=vec3(nx,ny,sz);
	z2*=rot2;
	
	sx=z2.x;
	sy=z2.y;
	sz=z2.z;
	nx=sx*sr12+sy*sr12;
	sy=-sx*sr12+sy*sr12; 
	sx=nx;
	nx=sx*sr23+sz*sr13;
	sz=-sx*sr13+sz*sr23;				//some things can be cleaned up
	sx=nx;
	
							//  Here is the end of the section
							//  the above code can be used with other fractal types

	if ((manPower==2.)) {
		float sx2=sx*sx;
		float sy2=sy*sy;				// will be switching code around later       
		float sz2=sz*sz;
		nx=sx2-sy2-sz2;   //
		r1=2.*(sx)/sqrt(sy2+sz2);  //  you can switch ny and nz... I'll add a boolean later
		float nz=r1*(sy2-sz2);
 		ny=r1*2.*sy*sz;   //
		z= vec3(nx,ny,nz);  
	} else {
		rxyz=pow((sx*sx+sy*sy+sz*sz),abs(manPower)/2.);
		r1=sqrt(sy*sy+sz*sz);
		victor=atan(r1,sx)*manPower;
		float phi=atan(sz,sy)*manPower;
		sx=rxyz*cos(victor);
		r1=rxyz*sin(victor);
		sz=r1*cos(phi);
		sy=r1*sin(phi);
		z= vec3(sx,sy,sz);
	}
}


bool inside(vec3 pos) {
	
	//vec3 z=pos;
		float ryz=sqrt(pos.y*pos.y+pos.z*pos.z);
		float theta=(Spokes)*atan(pos.z,pos.y)/3.;
		vec3 z=vec3(pos.x,ryz*cos(theta),ryz*sin(theta));
	 
	float r;
	int i=0;
	r=length(z);
	float d=magPower+magSoft;
	while(r<Bailout && (i<Iterations)) {
		
			magpart(z,d);
			z*=rot1;
			
		if (i>ColorIterationsMin && i<ColorIterations) {
			orbitTrap=min(orbitTrap,vec4(z.x*z.x,z.y*z.y,z.z*z.z,dot(z,z)));

		}
		if (Julia) {
			z+=JuliaC;	
		} else {
			if (pos.x>0.) {                          // can I do this???
				z+=vec3(pos.x*.5,0.,0.);
			} else {
				z+=vec3(pos.x,0.,0.);
			}
		}	
		r=length(z);
		z*=rot3;
		i++;
	}
	return (r<Bailout);
}




#preset default
FOV = 0.62536
Eye = 2.02493,-0.000258894,0.00365772
Target = -0.498369,0.00627817,-0.0886997
Up = -0.0024414,-0.999997,0
EquiRectangular = false
Gamma = 2.5
ToneMapping = 3
Exposure = 1.34694
Brightness = 1
Contrast = 0.9901
Saturation = 1
Near = 0.48984
Far = 2.71704
NormalScale = 0.00024
AOScale = 0.00291
Glow = 0.34167
AOStrength = 0.86047
Samples = 254
Stratify = true
DebugInside = false
SampleNeighbors = true
Specular = 1
SpecularExp = 5.455
SpotLight = 0.027451,0.196078,0.658824,0.89706
SpotLightDir = 0.80246,0.78126
CamLight = 1,1,1,1.44928
CamLightMin = 0.54217
Fog = 1.84
ShowDepth = false
DebugNormals = false
BaseColor = 0.388235,0.388235,0.388235
OrbitStrength = 0.78
X = 0.329412,0.321569,0.364706,0.08738
Y = 0.262745,0.333333,0.282353,0.02912
Z = 0.054902,0.0431373,0.439216,0.1068
R = 0.317647,0.670588,0.690196,0.64706
BackgroundColor = 0.0196078,0.0117647,0.192157
GradientBackground = 1.73915
CycleColors = false
Cycles = 23.1868
Iterations = 6
ColorIterations = 4
Bailout = 12
Julia = false
JuliaC = 0,0,0
ColorIterationsMin = 2
#endpreset

  Should clean up that code.   Anyways- couldn't get animations to work, for whatever reason, so I had to use ChaosPro for the following (also, am going to post a "pseudopalette" type idea I have for Framentarium, sometime soon- basically use a jpeg of a palette and Texture2D calls... maybe can do the palette type stuff that I love- already messed with textures earlier- nothing too great as of yet).

  Here are the gears turning in my brain, for those of you with a sense of humor, you might notice the slippage from my youtube splice (without deleting the extra frames).  It's Rotation2, vector (1,1,1) in the above formula:
<a href="http://www.youtube.com/v/Ih4k5KZAVxg&rel=1&fs=1&hd=1" target="_blank">http://www.youtube.com/v/Ih4k5KZAVxg&rel=1&fs=1&hd=1</a>

[richardrosenman]

Your animation looks awesome, Benesi - nice job!

-Rich