# Author: Jiwon Hahn
# ECE 207 Spring '02
# Generate an image using ray tracing for a simple scene


from math import *
from array import *	#to read in the file values as 'unsigned char'


# define values
r0 = 200     # radius of the metal sphere
z0 = 200     # the center coordinate of the sphere 
z = -50     # distance from the image to the picture
W=200  # width of the image 
H=150  # height of the image

# funcion: normalize V to unit vector
def norm(V):
	s=sqrt(V[0]*V[0]+V[1]*V[1]+V[2]*V[2])
        V[0], V[1], V[2] = V[0]/s,V[1]/s,V[2]/s
        return V

# initialize arrays with type 'unsigned char'='B'
a1,a2,a3=array('B',[]),array('B',[]),array('B',[])
	
# open r,g,b files to read.
f1=open('cat.b','r')
f2=open('cat.g','r')
f3=open('cat.r','r')

# read the files and close.
a1.fromfile(f1,307200) #640x480
a2.fromfile(f2,307200)
a3.fromfile(f3,307200)
f1.close()
f2.close()
f3.close()

blines,glines,rlines=a1.tolist(),a2.tolist(),a3.tolist()
print len(glines) # test

### store the r,g,b values in 2D arrays ###
picr,picg,picb=[],[],[]
tr,tg,tb=[],[],[] #temp storage for 1D array elements
k=0

for i in range(480):    #initialization of arrays
    picr.append(0)
    picg.append(0)
    picb.append(0)

for pi in range(480):   #creating 2D arrays to store picture values 
    for pj in range(640):
        tr.append(rlines[k])
        tg.append(glines[k])
        tb.append(blines[k])
        k=k+1
    picr[pi],picg[pi],picb[pi]=tr,tg,tb
    tr,tg,tb=[],[],[]

# create & initialize 2D arrays to store image values
imgr,imgg,imgb=[],[],[]

for i in range(H):
    imgr.append(0)
    imgg.append(0)        
    imgb.append(0)

for i in range(H):
    for j in range(W):
        tr.append(0)
	tg.append(0)
	tb.append(0)
    imgr[i],imgg[i],imgb[i]=tr,tg,tb
    tr,tg,tb=[],[],[]
	 

### Main Algorithm ###
for i in range(H):
    for j in range(W):
        x=W/2-j # set x,y,z<0   
        y=i-H/2
        s=z*z*z0*z0-(x*x+y*y+z*z)*(z0*z0-r0*r0)	#inside the square root


	########## 1.Ray doesn't meet the sphere (blue) ##########
	if(s<=0):	
            imgr[i][j],imgg[i][j],imgb[i][j]=0,0,255
 
	else:
	    t=(z*z0-sqrt(s))/(x*x+y*y+z*z)
	    X,Y,Z = x*t,y*t,z*t   #intersecting point on sphere

	    I=norm([x,y,z])	  #sphere to image vectoR 
            N=norm([X,Y,Z-z0])    #normal vector of the sphere
	    IN=I[0]*N[0]+I[1]*N[1]+I[2]*N[2]
	    P=norm([2*IN*N[0]-I[0],2*IN*N[1]-I[1],2*IN*N[2]-I[2]]) #sphere to picture vector

	    t=Z/P[2]
            px,py = X-P[0]*t, Y-P[1]*t 	#(px,py)=picture coordinate. 


	    ########## 2.The ray goes to the opposite direction(ambient) ##########
	    if(P[2]>=0):
               		imgr[i][j],imgg[i][j],imgb[i][j]=255,255,255      

       
	    ########## 3. The ray traced from the image hits the picture (cat) ##########
            elif((-320<=px)and(px<320)and(-240<=py)and(py<240)):
                	pi,pj=int(240-py),int(320+px)
                 	imgr[i][j],imgg[i][j],imgb[i][j]= \
				 picr[pi][pj],picg[pi][pj],picb[pi][pj]


	    ########## 4. Out of picture range (ambient) ##########
            else:   	
               		imgr[i][j],imgg[i][j],imgb[i][j]=255,255,255                                                          


# Convert 2D arrays to 1D arrays.
rLines,gLines,bLines=[],[],[]
for i in range(H):
    rLines=rLines+imgr[i]
    gLines=gLines+imgg[i]
    bLines=bLines+imgb[i]

# Convert from int to (unsigned)char
rLines=map(chr,rLines)
gLines=map(chr,gLines)
bLines=map(chr,bLines)

# Write to r,b,g files.
f1=open('out.r','w')
f2=open('out.g','w')
f3=open('out.b','w')
f1.writelines(rLines)
f2.writelines(gLines)
f3.writelines(bLines)

f1.close()
f2.close()
f3.close()