// 'Ponder' by Mike Andrews, 3rd September 2006
// Figure not included ...

// +h600 +w971 +a0.2 +am2 +r2

// --- code starts ---

// --- MegaPov 1.2.1 only! ---
#version unofficial megapov 1.21;

#declare _Phi=(sqrt(5)-1)/2;
#declare Radiosity = on;
#declare ModifySamples = on;
#declare Figure = off;

global_settings{
  exposure 1 exposure_gain 1
  tone_mapping {
    function { ln(x+1) }
    inverse function { exp(x)-1 }
  }
  #if(Radiosity)
    radiosity {
      randomize on
      #declare maxS = 100;
      #if (ModifySamples)
        #declare Spread = sin(radians(60)); // normally 90
        #declare S = array[maxS];
        #declare RS = seed(87576);
        samples {
          maxS,
          #declare C = 0; #while (C < maxS)
            #declare F = pow((C+0.5)/maxS,1/3); // normally 1/2
            #declare S[C] = vrotate(x, _Phi*C*360*z)*F*Spread;
            #declare S[C] = <S[C].x, S[C].y, sqrt(1-S[C].x*S[C].x-S[C].y*S[C].y)>;
          #declare C = C + 1; #end
          #declare C = 0; #while (C < maxS)
            #declare N = rand(RS);
            #declare T = S[C];
            #declare S[C] = S[N];
            #declare S[N] = T;
          #declare C = C + 1; #end
          #declare C = 0; #while (C < maxS)
            <S[C].x, S[C].y, S[C].z> #if (C < maxS-1) , #end
          #declare C = C + 1; #end
        }
      #end
      pretrace_start 8/image_width
      pretrace_end   2/image_width
      count maxS
      nearest_count 20
      error_bound 0.05
      recursion_limit 3
      low_error_factor .02
      gray_threshold 0.0
      minimum_reuse 0.001
      brightness 1
      adc_bailout 0.001/2
      normal on
    }
  #end
  assumed_gamma 1
  max_trace_level 20
}

#macro vPow(_V,_P) <pow(_V.x,_P),pow(_V.y,_P),pow(_V.z,_P)> #end

#declare BC = <4/7,5.5/7,1>;

sky_sphere { 
  pigment { 
    planar 
    colour_map {
      [0  rgb 2*vPow(BC,4)]
      [.9 rgb 2]
      [1  rgb 2*(1-vPow(BC,4))]
    } 
  } 
}

#declare camPos = 9*<-1,3,-5>;

#declare Phi=(sqrt(5)+1)/2;
#declare rPhi=2-Phi;

#declare Cen = (Phi-1)/(2*Phi-1);

camera{
  up y
  right x*image_width/image_height
  location camPos
  look_at  -y*rPhi*20
  angle 40
}

#declare Ill = box {
  <0,-_Phi+0.01,1.01>,<0.02,-rPhi-0.01,Phi-0.01>
  translate <-Cen,0,-Phi*(1-Cen)>
}

union {
  #declare C = 0; #while (C < 10)
    object { 
      Ill scale pow(_Phi,C) rotate -90*C*y 
      pigment { rgb <.8,.5,.2> }
      finish { diffuse 0 ambient 10 }
    }
    #debug concat("object ", str(C,0,0),"\n")
  #declare C = C + 1; #end
  translate <Cen,0,Phi*(1-Cen)>-<0.5,0,Phi/2>
  rotate 90*y
  scale 20
}

#declare Cut = difference {
  box { -<0,_Phi,0>,<1,-rPhi,Phi> }
  box { <-0.01,-_Phi+0.01,1.01>,<0.01,-rPhi-0.01,Phi-0.01> }
  translate <-Cen,0,-Phi*(1-Cen)>
}

#if (Figure) #include "vp_POV_geom.inc" #end

union {
  plane { y, 1-Phi }
  #declare C = 0; #while (C < 10)
    object { Cut scale pow(_Phi,C) rotate -90*C*y }
    sphere { 0, 0.1 translate -<Cen,rPhi,Phi*(1-Cen)>+<0.5,0.1,0.5> scale pow(_Phi,C) rotate -90*C*y }
    union {
      #declare P = vnormalize(-<Cen,rPhi,Phi*(1-Cen)>)*(vlength(<Cen,rPhi,Phi*(1-Cen)>)-0.04);
      #declare nD = 13;
      #declare D = 0.25; #while (D < nD)
        #declare S1 = pow(_Phi, D/nD);
        #declare P1 = vrotate(P*S1,-90*D*y/nD);
        sphere { 
          P1, 0.04*S1 
          pigment { rgb <.8,.5,.2> }
          finish { diffuse 0.3 reflection { 0.2,0.6 metallic } metallic ambient 0 }
        }
        cylinder { <P1.x,-rPhi-0.1,P1.z>,P1,0.04*_Phi*S1 }
      #declare D = D + 1; #end
      scale pow(_Phi,C) rotate -90*C*y
    }
    #debug concat("object ", str(C,0,0),"\n")
  #declare C = C + 1; #end
  #if (Figure)
    object {
      vp_
      hollow
      double_illuminate
      rotate -30*y
      scale 0.2/0.53
      translate -<Cen,rPhi,Phi*(1-Cen)>+<0.5,0,0.5>
    }
  #end
  translate <Cen,0,Phi*(1-Cen)>-<0.5,0,Phi/2>
  rotate 90*y
  scale 20
  pigment { rgb 1 }
  finish { diffuse 0.25 ambient 0 }
}

// --- code ends ---