/*
Pov3 Parametric equation - Moebius strip.
By Tore Nordstrand  nfytn@alf.uib.no
http://www.uib.no/People/nfytn/mathgal.htm
*/
#include "colors.inc"
#version 3.0
global_settings { assumed_gamma 1.0 }
#declare net   =  1           //0 = triangles, 1 = ribs.
#declare uiter = 28            //number of steps (lengthwise)
#declare viter = 14            //and across. Edit these for better
                               //resolution

#declare vmin  = -0.3         //endpoints in the (width)
#declare vmax  =  0.3         //parametric rectangle
#declare umin  =  0   
#declare umax  =  2*pi
#declare iu    = (umax-umin)/uiter    //increment per step
#declare iv    = (vmax-vmin)/viter 
//Simple textures - edit these for better results
#declare simple1 =  texture{
                      pigment {color Red}
                      finish {ambient 0.2 specular 0.7 roughness 0.05 }
                           }
#declare simple2 =  texture{
                      pigment {color Green}
                      finish {ambient 0.2 specular 0.7 roughness 0.05 }
                           }

//Lights, camera, action
camera
{ location  < 15, 9,-15>
  direction  1*z
  look_at   < 0, 0, 0>
} 

light_source { <0, 0,-19> color rgb 1 }
light_source { <23, 20, -15> color rgb 1 }
//object {box{-30,30 inverse}pigment{color rgb 0.8}}



//STRIP STARTS HERE

union{
  #declare uu = umin    
  #while (uu<=umax)                        //outer loop
    #declare vv = vmin
    #while (vv<=vmax)                      //inner loop

//POINTS
      #declare x1=cos(uu)+vv*cos(uu/2)*cos(uu) //Calculate the
      #declare y1=sin(uu)+vv*cos(uu/2)*sin(uu) //four points in
      #declare z1=vv*sin(uu/2)                 //the rectangle defined
      #declare uu=uu+iu                        //by uu,uu+iu,vv,vv+iv
      #declare x2=cos(uu)+vv*cos(uu/2)*cos(uu) //Uses uu and vv because
      #declare y2=sin(uu)+vv*cos(uu/2)*sin(uu) //u and v are built in
      #declare z2=vv*sin(uu/2)                 //constants.
      #declare vv=vv+iv  
      #declare x3=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare y3=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare z3=vv*sin(uu/2)
      #declare uu=uu-iu
      #declare x4=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare y4=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare z4=vv*sin(uu/2)
      #declare vv=vv-iv
#if (net)
#else     //do this only for trangles
//NORMALS   used for smoothing triangles
      #declare vv=vv-iv                         //Calculate surface normals
      #declare nx1=cos(uu)+vv*cos(uu/2)*cos(uu) //to the nine rectangles in
      #declare ny1=sin(uu)+vv*cos(uu/2)*sin(uu) //the immediate neighbourhood 
      #declare nz1=vv*sin(uu/2)
      #declare uu=uu+iu 
      #declare nx2=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny2=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz2=vv*sin(uu/2)
      #declare uu=uu+iu 
      #declare vv=vv+iv
      #declare nx3=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny3=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz3=vv*sin(uu/2)
      #declare vv=vv+iv
      #declare nx4=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny4=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz4=vv*sin(uu/2)
      #declare vv=vv+iv 
      #declare uu=uu-iu
      #declare nx5=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny5=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz5=vv*sin(uu/2)
      #declare uu=uu-iu
      #declare nx6=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny6=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz6=vv*sin(uu/2)
      #declare uu=uu-iu
      #declare vv=vv-iv
      #declare nx7=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny7=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz7=vv*sin(uu/2)
      #declare vv=vv-iv
      #declare nx8=cos(uu)+vv*cos(uu/2)*cos(uu)
      #declare ny8=sin(uu)+vv*cos(uu/2)*sin(uu)
      #declare nz8=vv*sin(uu/2)
      #declare uu=uu+iu
      //here is where normals are calculated.
      #declare n1=vcross(<x1-nx1,y1-ny1,z1-nz1>,<x1-nx8,y1-ny8,z1-nz8>)
      #declare n2=vcross(<x1- x2,y1- y2,z1- z2>,<x1-nx1,y1-ny1,z1-nz1>)   
      #declare n3=vcross(<x1- x4,y1- y4,z1- z4>,<x1- x2,y1- y2,z1- z2>)   
      #declare n4=vcross(<x1-nx8,y1-ny8,z1-nz8>,<x1- x4,y1- y4,z1- z4>)
      #declare normal1=n1+n2+n3+n4
      #declare n1=vcross(<x2-nx2,y2-ny2,z2-nz2>,<x2- x1,y2- y1,z2- z1>)
      #declare n2=vcross(<x2-nx3,y2-ny3,z2-nz3>,<x2-nx2,y2-ny2,z2-nz2>)   
      #declare n3=vcross(<x2- x3,y2- y3,z2- z3>,<x2-nx3,y2-ny3,z2-nz3>)   
      #declare n4=vcross(<x2- x1,y2- y1,z2- z1>,<x2- x3,y2- y3,z2- z3>)
      #declare normal2=n1+n2+n3+n4
      #declare n1=vcross(<x3- x2,y3- y2,z3- z2>,<x3- x4,y3- y4,z3- z4>)
      #declare n2=vcross(<x3-nx4,y3-ny4,z3-nz4>,<x3- x2,y3- y2,z3- z2>)   
      #declare n3=vcross(<x3-nx5,y3-ny5,z3-nz5>,<x3-nx4,y3-ny4,z3-nz4>)   
      #declare n4=vcross(<x3- x4,y3- y4,z3- z4>,<x3-nx5,y3-ny5,z3-nz5>)
      #declare normal3=n1+n2+n3+n4
      #declare n1=vcross(<x4- x1,y4- y1,z4- z1>,<x4-nx7,y4-ny7,z4-nz7>)
      #declare n2=vcross(<x4- x3,y4- y3,z4- z3>,<x4- x1,y4- y1,z4- z1>)   
      #declare n3=vcross(<x4-nx6,y4-ny6,z4-nz6>,<x4- x3,y4- y3,z4- z3>)   
      #declare n4=vcross(<x4-nx7,y4-ny7,z4-nz7>,<x4-nx6,y4-ny6,z4-nz6>)
      #declare normal4=n1+n2+n3+n4
#end //if

#if (net)
//CYLINDERS
      object{cylinder{<x1,y1,z1>,<x2,y2,z2>,0.013}
             texture{simple1}
             }
      object{cylinder{<x3,y3,z3>,<x4,y4,z4>,0.013}
             texture{simple1}
             }

      object{cylinder{<x1,y1,z1>,<x4,y4,z4>,0.013}
             texture{simple2}
             }

#else
//TRIANGLES
      object{smooth_triangle{<x1,y1,z1>, normal1,
                             <x2,y2,z2>, normal2,
                             <x3,y3,z3>, normal3}
               texture{simple1}
             }
      object{smooth_triangle{<x1,y1,z1>, normal1,  
                             <x3,y3,z3>, normal3,   
                             <x4,y4,z4>, normal4}
               texture{simple1}
              }

#end //object if

      #declare vv = vv+iv
   #end                 //inner loop
   #declare uu = uu+iu
 #end                   //outer loop
rotate <90,0,0>     //This is just changing the
rotate <0,90,0>     //coordinate system.
scale 9
}