Homework 0: 2D光线追踪 vs JavaScript

同样场景,Taichi版比p5.js版大概快了20倍左右(!)
现在还只支持线段与射线相交,而且衰减函数也是拍脑袋想的,而且可能是因为我写的有误差再加上一些微妙的差别,两个版本的结果看起来有些不同,不过大体上我已经知道了我想知道的这两者之间的数量级的差距 XD

JavaScript版:
https://editor.p5js.org/nitroglycerine/sketches/yDq3_-hY2

Taichi版:

# -*- coding: utf-8 -*-
import taichi as ti

# Javascript:         0.1  FPS
# i7-6820HQ     CPU:  2.70 FPS
# Quadro M1000M CUDA: 2.03 FPS
ti.init(arch=ti.cpu)

# 把场景写进kernel导致编译时间太长,故关掉高级优化
ti.core.toggle_advanced_optimization(False)

W = 720; H = 480;
MAGNIFY_DISP = 4
color_buffer     = ti.Vector(3, dt=ti.f32, shape=(W, H));
gui = ti.GUI('Haha', (W, H))
# Taichi's coordinate system is different from Processing's
# 
#     Taichi                   Processing
#                                  
# ^ (0, 480)                +------------> (720,0
# |                         |
# |                         |
# +----------> (720, 0)     V (0, 480)

@ti.func
def Perp(v):
  return ti.Vector([-v[1], v[0], 0.0])

# 返回:[t, tact]
@ti.func
def IntersectLineWithRay(l, r):
  p1 = l[0]; p2 = l[1]; o = r[0]; d = ti.normalized(r[1]);
  p1p2 = p2 - p1
  n = ti.normalized(Perp(p1p2))
  op1 = p1 - o;
  dist = ti.dot(n, op1);
  ddn  = ti.dot(n, d);
  t    = dist / ddn;
  tact = ti.Vector([0.0, 0.0, 0.0]);
  if t > 0:
    tact = o + t * d
    t0 = ti.dot((tact - p1), p1p2)
    t1 = ti.dot((p2 - tact), p1p2)
    completion = t0 / (t0+t1)
    if completion >= 0 and completion <= 1:
      pass
    else:
      t = -1
  return t, tact 

@ti.kernel
def Clear():
  for i, j in color_buffer:
    color_buffer[i, j] = ti.Vector([0,0,0]);

@ti.func
def Diffuse(l, d):
  n = Perp(l[1]-l[0]);
  if (ti.dot(n, d) > 0): n = n*(-1)
  theta = ti.random(ti.f32) * 3.14159*2
  ret = ti.Vector([ti.cos(theta), ti.sin(theta), 0.0]);
  for x in range(0, 4):
    theta = ti.random(ti.f32) * 3.14159*2
    ret = ti.Vector([ti.cos(theta), ti.sin(theta), 0.0]);
    if ti.dot(n, ret) < 0:
      break
  return ret

NSAMP = 4;

x0 = 0.2*W; x1 = 0.25*W; x2 = 0.3*W; y1 = 0.7*H; y2 = 0.55*H; y0 = 0.85*H;
g1 = [ti.Vector([x1,y1,0]), ti.Vector([x2,y1,0])]
g2 = [ti.Vector([x2,y1,0]), ti.Vector([x2,y2,0])]
g3 = [ti.Vector([x0,y2,0]), ti.Vector([x2,y2,0])]
g4 = [ti.Vector([x0,y2,0]), ti.Vector([x0,y0,0])]
g5 = [ti.Vector([x0,y0,0]), ti.Vector([x2,y0,0])]

#####
#
# ###
#   #
#####
@ti.func
def IntersectG(o, d, t_, this_intensity_, done_):
  d_out = d;
  t = t_;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(g1, [o, d]);
  if t1 > 0 and t>t1:
    t = t1;
    this_intensity = ti.Vector([0.5, 0.0, 0.0]);
    done = True
  t1, _ = IntersectLineWithRay(g2, [o, d]);
  if t1 > 0 and t>t1: t = t1; d_out = Diffuse(g2, d); done = False;
  t1, _ = IntersectLineWithRay(g3, [o, d]);
  if t1 > 0 and t>t1: t = t1; this_intensity = ti.Vector([0.5, 0.0, 0.0]); done = True
  t1, _ = IntersectLineWithRay(g4, [o, d]);
  if t1 > 0 and t>t1: t = t1; d_out = Diffuse(g4, d); done = False
  t1, _ = IntersectLineWithRay(g5, [o, d]);
  if t1 > 0 and t>t1: t = t1; this_intensity = ti.Vector([0.5, 0.0, 0.0]); done = True
  return t, d_out, this_intensity, done

  #
 # #
#   #
#####
#   #
x0 = 0.325*W; x1=0.375*W; x2=0.425*W; y1=0.7*H; y2=0.55*H; y0=0.85*H; y3=0.65*H;
a1 = [ti.Vector([x1,y0,0]), ti.Vector([x2,y1,0])]
a2 = [ti.Vector([x1,y0,0]), ti.Vector([x0,y1,0])]
a3 = [ti.Vector([x0,y3,0]), ti.Vector([x2,y3,0])]
a4 = [ti.Vector([x0,y1,0]), ti.Vector([x0,y2,0])]
a5 = [ti.Vector([x2,y1,0]), ti.Vector([x2,y2,0])]

@ti.func
def IntersectA(o, d, t_, this_intensity_, done_):
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(a1, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(a1, d); done = False;
  t1, _ = IntersectLineWithRay(a2, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(a2, d); done = False;
  t1, _ = IntersectLineWithRay(a3, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.0, 0.8, 0.0]); done = True
  t1, _ = IntersectLineWithRay(a4, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(a4, d); done = False;
  t1, _ = IntersectLineWithRay(a5, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(a5, d); done = False;
  return t, d_out, this_intensity, done;
  
#   #
## ##
# # #
#   #
#   #
x0=0.45*W; x1=0.5*W; x2=0.55*W; y1=0.7*H; y2=0.55*H; y0=0.85*H;
m1 = [ti.Vector([x0,y0,0]), ti.Vector([x0,y2,0])]
m2 = [ti.Vector([x0,y0,0]), ti.Vector([x1,y2,0])]
m3 = [ti.Vector([x2,y0,0]), ti.Vector([x1,y2,0])]
m4 = [ti.Vector([x2,y0,0]), ti.Vector([x2,y2,0])]

@ti.func
def IntersectM(o, d, t_, this_intensity_, done_):
  
  d_out = d;
  t = t_;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(m1, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(m1, d); done = False;
  t1, _ = IntersectLineWithRay(m2, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.4, 0.4, 0.0]); done = True
  t1, _ = IntersectLineWithRay(m3, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.4, 0.4, 0.0]); done = True
  t1, _ = IntersectLineWithRay(m4, [o, d]);
  if t1 > 0 and t > t1: t = t1;  d_out = Diffuse(m4, d); done = False;
  return t, d_out, this_intensity, done;

#####
#
#####
#
#####
x0=0.575*W; x1=0.625*W; x2=0.675*W; y1=0.7*H; y2=0.55*H; y0=0.85*H;
e1 = [ti.Vector([x0,y0,0]), ti.Vector([x0,y2,0])]
e2 = [ti.Vector([x0,y0,0]), ti.Vector([x2,y0,0])]
e3 = [ti.Vector([x0,y2,0]), ti.Vector([x2,y2,0])]
e4 = [ti.Vector([x0,y1,0]), ti.Vector([x1,y1,0])]

@ti.func
def IntersectE(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(e1, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(e1, d); done = False;
  t1, _ = IntersectLineWithRay(e2, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0, 0.5, 0.5]); done = True
  t1, _ = IntersectLineWithRay(e3, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0, 0.5, 0.5]); done = True
  t1, _ = IntersectLineWithRay(e4, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0, 0.5, 0.5]); done = True
  return t, d_out, this_intensity, done;

#####
#
#####
    #
#####
x0=0.70*W; x1=0.75*W; x2=0.80*W; y1=0.7*H; y2=0.55*H; y0=0.85*H;
s1 = [ti.Vector([x0,y0,0]), ti.Vector([x2,y0,0])]
s2 = [ti.Vector([x0,y0,0]), ti.Vector([x0,y1,0])]
s3 = [ti.Vector([x0,y1,0]), ti.Vector([x2,y1,0])]
s4 = [ti.Vector([x2,y1,0]), ti.Vector([x2,y2,0])]
s5 = [ti.Vector([x0,y2,0]), ti.Vector([x2,y2,0])]
@ti.func
def IntersectS(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_
  done = done_
  t1, _ = IntersectLineWithRay(s1, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.5, 0.2, 0.2]); done = True
  t1, _ = IntersectLineWithRay(s2, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(s2, d); done = False;
  t1, _ = IntersectLineWithRay(s3, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.5, 0.2, 0.2]); done = True
  t1, _ = IntersectLineWithRay(s4, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(s4, d); done = False;
  t1, _ = IntersectLineWithRay(s5, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.5, 0.2, 0.2]); done = True
  return t, d_out, this_intensity, done;

#####
    #
#####
#
#####
x0=0.25*W; x1=0.3*W; x2=0.35*W; y1=0.3*H; y2=0.15*H; y0=0.45*H;
v21 = [ti.Vector([x0,y0,0]), ti.Vector([x2,y0,0])]
v22 = [ti.Vector([x2,y0,0]), ti.Vector([x2,y1,0])]
v23 = [ti.Vector([x0,y1,0]), ti.Vector([x2,y1,0])]
v24 = [ti.Vector([x0,y1,0]), ti.Vector([x0,y2,0])]
v25 = [ti.Vector([x0,y2,0]), ti.Vector([x2,y2,0])]
@ti.func
def Intersect2(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(v21, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v21, d); done = False;
  t1, _ = IntersectLineWithRay(v22, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v22, d); done = False;
  t1, _ = IntersectLineWithRay(v23, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v23, d); done = False;
  t1, _ = IntersectLineWithRay(v24, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v24, d); done = False;
  t1, _ = IntersectLineWithRay(v25, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v25, d); done = False;
  return t, d_out, this_intensity, done;
  
#####
#   #
#   #
#   #
#####
x0=0.45*W; x1=0.5*W; x2=0.55*W; y1=0.3*H; y2=0.15*H; y0=0.45*H;
v01 = [ti.Vector([x0,y0,0]), ti.Vector([x2,y0,0])]
v02 = [ti.Vector([x2,y0,0]), ti.Vector([x2,y2,0])]
v03 = [ti.Vector([x2,y2,0]), ti.Vector([x0,y2,0])]
v04 = [ti.Vector([x0,y0,0]), ti.Vector([x0,y2,0])]
@ti.func
def Intersect0(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(v01, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.2, 0.2, 0.4]); done = True
  t1, _ = IntersectLineWithRay(v02, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v02, d); done = False;
  t1, _ = IntersectLineWithRay(v03, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.2, 0.2, 0.4]); done = True
  t1, _ = IntersectLineWithRay(v04, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v04, d); done = False;
  return t, d_out, this_intensity, done;
  
  #
 ##
# #
  #
#####
x0=0.65*W; x1=0.7*W; x2=0.75*W; y1=0.3*H; y2=0.15*H; y0=0.45*H;
v11 = [ti.Vector([x1,y0,0]), ti.Vector([x1-0.05*W,y0-0.075*H,0])]
v12 = [ti.Vector([x1,y0,0]), ti.Vector([x1,y2,0])]
v13 = [ti.Vector([x2,y2,0]), ti.Vector([x0,y2,0])]
@ti.func
def Intersect1(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(v11, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.4, 0.4, 0.8]); done = True
  t1, _ = IntersectLineWithRay(v12, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v12, d); done = False;
  t1, _ = IntersectLineWithRay(v13, [o, d]);
  if t1 > 0 and t > t1: t = t1; d_out = Diffuse(v13, d); done = False;
  return t, d_out, this_intensity, done;

##############
#            #
#            #

#            #
#            #
##############
vw1 = [ti.Vector([0.15*W,0.5*H,0]), ti.Vector([0.15*W,0.1*H,0])]
vw2 = [ti.Vector([0.85*W,0.5*H,0]), ti.Vector([0.85*W,0.1*H,0])]
vw3 = [ti.Vector([0.15*W,0.05*H,0]), ti.Vector([0.85*W,0.05*H,0])]
vw4 = [ti.Vector([0.15*W,0.95*H,0]), ti.Vector([0.85*W,0.95*H,0])]
vw5 = [ti.Vector([0.15*W,0.9*H,0]), ti.Vector([0.15*W,0.5*H,0])]
vw6 = [ti.Vector([0.85*W,0.9*H,0]), ti.Vector([0.85*W,0.5*H,0])]
@ti.func
def IntersectW(o, d, t_, this_intensity_, done_):
  
  t = t_;
  d_out = d;
  this_intensity = this_intensity_;
  done = done_
  t1, _ = IntersectLineWithRay(vw1, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.6, 0.6, 1.4]); done = True
  t1, _ = IntersectLineWithRay(vw2, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.6, 0.6, 1.4]); done = True
  t1, _ = IntersectLineWithRay(vw3, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.6, 0.6, 1.4]); done = True
  t1, _ = IntersectLineWithRay(vw4, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.4, 0.8, 0.8]); done = True
  t1, _ = IntersectLineWithRay(vw5, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.1, 0.2, 0.2]); done = True
  t1, _ = IntersectLineWithRay(vw6, [o, d]);
  if t1 > 0 and t > t1: t = t1; this_intensity = ti.Vector([0.1, 0.2, 0.2]); done = True
  
  return t, d_out, this_intensity, done;

@ti.kernel
def Fill():
  for i, j in color_buffer:
    intensity = ti.Vector([0.0, 0.0, 0.0])
    samp = 0
    while samp < NSAMP: # 每个像素采这么多次
      o = ti.Vector([i*1.0, j*1.0, 0.0]); # 必须是三维的不能是二维的?
      theta = ti.random(ti.f32) * 3.14159*2
      d = ti.Vector([ti.cos(theta), ti.sin(theta), 0.0])
      # obj = None 不能设为none?
      
      t_total = 0.0 # 光走过的距离
      
      level = 0
      while level < 3:
        d_orig = d;
        t = 1000000000.0
        
        # 求交点
        # this_intensity = ti.Vector([0,0,0]) # 错误:整型的不能赋浮点值
        this_intensity = ti.Vector([0.0, 0.0, 0.0]) 
        done = False
        t, d, this_intensity, done = IntersectG(o, d_orig, t, this_intensity, done); # G
        t, d, this_intensity, done = IntersectA(o, d_orig, t, this_intensity, done); # A
        t, d, this_intensity, done = IntersectM(o, d_orig, t, this_intensity, done); # M
        t, d, this_intensity, done = IntersectE(o, d_orig, t, this_intensity, done); # E
        t, d, this_intensity, done = IntersectS(o, d_orig, t, this_intensity, done); # S
        t, d, this_intensity, done = Intersect2(o, d_orig, t, this_intensity, done); # 2
        t, d, this_intensity, done = Intersect0(o, d_orig, t, this_intensity, done); # 0
        t, d, this_intensity, done = Intersect1(o, d_orig, t, this_intensity, done); # 1
        t, d, this_intensity, done = IntersectW(o, d_orig, t, this_intensity, done); # Walls
        
        # 总结一下求交点结果决定是再发射间接光还是退出
        if t != 1000000000.0:
          t_total += t
        damp = 1 / (1 + ti.exp(t_total/400))
        EPS = 0.001;
        o = o + d_orig * (t - EPS);
        if done: intensity += this_intensity * damp; break
        level += 1 # 一次递归结束了
      samp += 1; # 一次采样结束了
    color_buffer[i, j] += intensity * 1.0 / NSAMP

g_frame_count = 0;
while True:
  while gui.get_event(ti.GUI.PRESS):
    if gui.event.key == ti.GUI.ESCAPE:
      exit()
  #Clear();
  Fill();
  g_frame_count += 1;
  img = (color_buffer).to_numpy() * (1.0 / (g_frame_count)) * MAGNIFY_DISP;
  gui.set_image(img);
  gui.show()
7 Likes

酷!特意构造了一个GAMES201的场景实在是有心了

1 Like