Ray tracing Bug with Diffuse material

Question

I am trying to do ray tracing in python (following the tutorial given in Ray Tracing in a weekend). Basically I am shooting rays from eye and the recursively bounce around, each time they hit something they become weaker (actually this is reverse ray casting but you get the idea). The output I am getting is incorrect.

My output:

Expected output:

The shadows are messed up. What could be wrong?

My code:

import numpy as np
import sys
import random
from PIL import Image
from math import * 
from util import *

width  = 60
height = 60

samples = 20

#-----------------------------------------------------------------------

def reflected(vector, axis):
    return vector - axis * 2 * vector.dot(axis)


def RandomPointInSphere():
    p = None
    while True:
        p = Vector(random.uniform(0,1),random.uniform(0,1),random.uniform(0,1))*2 - Vector(1,1,1)
        if(p.dot(p) < 1):
            break
    return p

    
def GetNearestObject(objects, ray):
    nearest_obj = None
    min_hit = Intersection(None, INF, None, None)
    
    for obj in objects:
        hit = obj.intersect(ray)
        if(hit.distance < min_hit.distance):
            nearest_obj = obj
            min_hit = hit
    
    return min_hit
    

#-----------------------------------------------------------------------

def RayColor(objects, ray):
    # Part 1: Diffuse Material
    result = GetNearestObject(objects, ray)
    if(result.point != None):
        P = result.point
        N = result.normal
        E = RandomPointInSphere()
        target = P + N + E
        
        newRay = Ray(ray.origin, (target - ray.origin).normalize())
        return RayColor(objects, newRay)*0.5
    else:
        t = 0.5 * (ray.direction.y + 1.0);
        color = Vector(1.0, 1.0, 1.0)*(1.0 - t) + Vector(0.5, 0.7, 1.0)*t
        color.x = min(color.x,1.0)
        color.y = min(color.y,1.0)
        color.z = min(color.z,1.0)
        return Vector(1,1,1)
    
#-----------------------------------------------------------------------
def main():
    global bitmap
    
    eye    = Vector(0,0,1)
    ratio  = float(width) / height
    screen = Screen(-1, 1 / ratio, 1, -1 / ratio, 0)
    
    objects = []
    objects.append(Sphere(Vector(-0.2,0,-1),  0.7, Material(Vector(0.1,0,0),  Vector(0.7,0,0),    Vector(1,1,1), 100, 0.5)))
    
    objects.append(Sphere(Vector(0,-9000,0),  9000-0.7, Material(Vector(0.1,0.1,0.1),Vector(0.6,0.6,0.6),Vector(1,1,1), 100, 0.5)))
    
    light = Light(Vector(5,5,5), Material(Vector(1,1,1),Vector(1,1,1),Vector(1,1,1)))
    
    for frame in range(1):
    
        img    = Image.new(mode = "RGB", size=(width, height), color=Color.WHITE)
        bitmap = img.load() # create the pixel data

#--------------------------------------------------------------
#---                    YOUR CODE HERE                      ---
#--------------------------------------------------------------
        sys.setrecursionlimit(10000)
        
        #breakpoint()
        
        deltaX = (screen.right - screen.left)/(width-1)
        deltaY = (screen.top - screen.bottom)/(height-1)

        for y in range(height):
            for x in range(width):
                pixel     = Vector(screen.left+x*deltaX, screen.top-y*deltaY, screen.z)
                direction = (pixel - eye).normalize()
                pixelRay = Ray(eye, direction)                

                # Part 1: Diffuse Material
                color = Vector(0,0,0)
                for s in range(samples):
                    color += RayColor(objects, pixelRay)
                
                color *= 1.0/samples
                #color = Vector(sqrt(color.x), sqrt(color.y), sqrt(color.z))
                bitmap[x,y] = (int(color.x*255), int(color.y*255), int(color.z*255))

            print("progress: %d %%" % ((y+1)/height*100.0))
            
#--------------------------------------------------------------
#--------------------------------------------------------------
#--------------------------------------------------------------

        img.show()
        img.save("pic1.png")
        #img.save("images/fig" + f'{frame:06}' + ".png")
        #print("Saving ---> images/fig" + f'{frame:06}' + ".png")
        #img.close()

main()

util.py file:

from math import * 
from abc import ABC, abstractmethod

PI  = 3.14159265
INF = 9999999999

#-----------------------------------------------------------------------

class Color:
    BLACK = (0,0,0)
    RED   = (245,23,32)
    BLUE  = (46,139,192)
    GREEN = (24,165,88)
    PINK  = (250,38,160)
    GOLD  = (248,210,16)
    TEAL  = (43,124,133)
    WHITE = (255,255,255)
    CHROMA= (0,255,0)
    
#-----------------------------------------------------------------------

def dot(a, b):
    return a.x*b.x + a.y*b.y + a.z*b.z

def normalize(a):
    mag = a.magnitude()
    return Vector(a.x/mag,a.y/mag,a.z/mag)

class Vector:
    def __init__(self,x,y,z):
        self.x = x
        self.y = y
        self.z = z
        self.w = 1

    def toString(self):
        return "S:["+"{:.3f}".format(self.x)+", "+"{:.3f}".format(self.y)+", "+"{:.3f}".format(self.z)+"]"

    def dot(self, b):
        return self.x*b.x + self.y*b.y + self.z*b.z
    
    def cross(self, b):
        return Vector(self.y*b.z-self.z*b.y, self.z*b.x-self.x*b.z, self.x*b.y-self.y*b.x)
        
    def magnitude(self):
        return sqrt(self.x*self.x + self.y*self.y + self.z*self.z)

    def normalize(self):
        mag = self.magnitude()
        return Vector(self.x/mag, self.y/mag, self.z/mag)

    # Provide "overridden methods via the "__operation__" notation; allows you to do, for example, a+b, a-b, a*b
    def __add__(self, b):
        return Vector(self.x + b.x, self.y+b.y, self.z+b.z)

    def __sub__(self, b):
        return Vector(self.x-b.x, self.y-b.y, self.z-b.z)

    def __mul__(self, b):
        if type(b) == float or type(b) == int:
            return Vector(self.x*b, self.y*b, self.z*b)
        elif type(b) == Vector:
            return Vector(self.x*b.x, self.y*b.y, self.z*b.z)
        else:
            print(type(b))
            assert False

class Material:
    def __init__(self, ambient, diffuse, specular, shininess = 0.0, reflection = 0.0):
        self.ambient   = ambient    # 3D Vector
        self.diffuse   = diffuse    # 3D Vector
        self.specular  = specular   # 3D Vector
        self.shininess = shininess  #   0 - 100
        self.reflection= reflection # 0.0 - 1.0

class Light:
    def __init__(self, position, material):
        self.position = position
        self.material  = material

class Hittable(ABC):
    @abstractmethod
    def intersect(self, ray):
        pass

class Sphere(Hittable):
    def __init__(self,center, radius, material):
        self.center   = center
        self.radius   = radius
        self.material = material
        
    def intersect(self, ray):
        OC = ray.origin - self.center
        a = ray.direction.dot(ray.direction)
        b = 2.0 * ray.direction.dot(OC)
        c = OC.dot(OC) - self.radius*self.radius
        descriminant = b*b - 4*a*c
        if(descriminant > 0):
            t1 = (-b + sqrt(descriminant))/2
            t2 = (-b - sqrt(descriminant))/2
            if t1 > 0 and t2 > 0:
                dist = min(t1,t2)
                point = ray.pointAtParameter(dist)
                normal = (point - self.center).normalize()
                return Intersection(point, dist, normal, self)
        return Intersection(None, INF, None, None)
        
                
    def normal(self, b):
        return (b - self.center).normalize()


class Ray:
    def __init__(self, origin, direction):
        self.origin    = origin
        self.direction = direction

    def pointAtParameter(self, t):
        return self.origin + self.direction*t
    
    def toString(self):
        return "R->{ "+ self.origin.toString() + " -- " + self.direction.toString() +" }"
       
class Intersection:
    def __init__(self, point, distance, normal, object):
        self.point    = point
        self.distance = distance
        self.normal   = normal
        self.object   = object
        
        
class Screen:
    def __init__(self,  left, top, right, bottom, z):
        self.left   = left
        self.top    = top
        self.right  = right
        self.bottom = bottom
        self.z      = z
        
class Camera:
    def __init__(self):
        pass  
```

Answer 1

I see a couple of problems by scanning quickly over your code.

1. You are using ray.origin as new origin when scattering the ray off the surface, you should however use the hit point as origin of the new ray. This might also explain why you see a kind of "screen aligned" shadow.
2. The new ray direction is just the direction ray.origin to target.
3. Is result.point actually the center of the sphere or the point of intersection on the surface? I can't say for sure as I don't have the rest of your code.

There is a lot that is not right, or seemingly unintentional.

What you can do however to debug such things is to return a debug color like e.g. normal at point, depth (normalized to some maximum), scatter direction, etc.