<languageVersion: 1.0;>
/*
    PenroseTiles.pbk
    
    Author: Syntopia (Mikael Hvidtfeldt Christensen)
    This version: 18. september 2010
    
    Credits:
    --------
    This is a simple port of Tomasz Dobrowolski's (tomkh) "Procedural aperiodic fractals" EvalDraw script.
    It originates from this FractalForums thread: 
    http://www.fractalforums.com/new-theories-and-research/procedural-aperiodic-fractals/
    
    Copyright & License:
    --------------------
    You can use this code for whatever purpose, but be sure to acknowledge the original author of the algorithm (Tomasz Dobrowolski).
*/
kernel PenroseTiles
<   namespace : "Syntopia.Fractals";
    vendor : "Syntopia";
    version : 1;
    description : "Penrose Tile Generator.";
>

{
    output pixel4 dst;
    
    parameter int2 size
    <
        minValue:int2(100, 100); maxValue:int2(2000, 2000); defaultValue:int2(800, 600);
    >;
    
    parameter float2 center
    <
        minValue:float2(-2.0, -1.0); maxValue:float2(2.0, 1.0); defaultValue:float2(0.8, 0.4);
    >;
    
    parameter float zoom
    <
        minValue:0.0;maxValue:20.0;defaultValue:1.4;
    >;
    
    
    parameter int iterations
    <
        minValue:1; maxValue:50; defaultValue:10;
    >;
    
    region generated()
    {
        return region(float4(0, 0, size.x, size.y));
    }
    
    dependent float pi,sc,d1,d2,a1;    
    dependent float2x2 m1,m2,m3,m4,m5;
    dependent float2 p1;
       
    void evaluateDependents()
	{
        pi = 3.1415 ; // must be built-in somewhere
        sc = 2.0/(sqrt(5.0)-1.0);  // inflation scale

        // transformations constants
        d1 = tan(54.0*pi/180.0);
        d2 = tan(18.0*pi/180.0);
        a1 = .5/cos(36.0*pi/180.0);
        float a2 = (1.0+a1)*.5;
        float a3 = tan(36.0*pi/180.0)*a2; 
        float cos1 = cos(144.0*pi/180.0)*sc;
        float  sin1 = sin(144.0*pi/180.0)*sc;
        float cos2 = cos(108.0*pi/180.0)*sc;
        float sin2 = sin(108.0*pi/180.0)*sc;
        
        m1 = float2x2(-sc,0.0,  0.0,sc); 
        p1 = float2(-a2,-a3);
        m2= float2x2(cos1,-sin1,sin1,cos1); 
        m3= float2x2(cos1,sin1,-sin1,cos1); 
        m4= float2x2(-cos2,sin2,sin2,cos2); 
        m5= float2x2(cos2,sin2,-sin2,cos2); 
    }
    
    pixel4 getValue(float2 z) {
        int triangleType = 0; 
        for(int k=0; k<iterations; k++) {
           if (triangleType == 0) {
              if (1.0 - d1*z.y - z.x > 0.0) {         
                 z *= m1; z.x += sc;
              } else if (1.0 - d2*z.y - z.x > 0.0) {
                 z += p1; z *= m2; triangleType = 1;          
              } else {      
                 z.x-=(1.0+a1); z*= m3;
              }      
           } else {
              if (d1*z.y - z.x > 0.0) { 
                z*=m4; triangleType = 0; 
              } else { 
                 z.x -= a1; z*= m5;         
              }
           }   
        }
        return (triangleType == 0) ? pixel4(218.0/256.0,213.0/256.0,118.0/256.0,1.0) :  pixel4(74.0/256.0,101.0/256.0,166.0/256.0, 1.0)  ;
    }
    
    // Convert from screen space coords to model coords.
    float2 toModelSpace(float2 p) {
        // TODO: Some of this stuff could be precalculated as dependents...
        float x0 = center.x;
        float y0 = center.y;
        float xSpan = 1.0/exp(zoom);
        float x1 = x0 - xSpan*0.5;	
        float x2 = x0 + xSpan*0.5;	
        float ySpan = xSpan * (float(size.y) / float(size.x));
        float y1 = y0 - ySpan*0.5;
        float y2 = y0 + ySpan*0.5;
        return float2(x1, y1) + p * float2(xSpan / float(size.x), ySpan / float(size.y));
    }
    
    void evaluatePixel()
    {
        // Simple anti-alias. aaSteps=1 corresponds to 3x3 sampling.
        int aaSteps = 1;
        float stepSize = 1.0/(3.0*float(aaSteps));
        dst = pixel4(0.0,0.0,0.0,0.0);
        for (int x = -aaSteps; x<=aaSteps; x++) {
            for (int y = -aaSteps; y<=aaSteps; y++) {
                dst += getValue(toModelSpace(outCoord()+float2(float(x)*stepSize,float(y)*stepSize)));
            }
        }
        dst /= float((aaSteps*2+1)*(aaSteps*2+1));
    }
}

PenroseTiles.pbk

C