Experiments In Ray Tracing 2 (materials, shaders)

[Hades]

Here is part 2!

There is a PB version of the intersection Procedure, too. (Just swap the '1' in the name between the two)

Edit: fixed a few problems, should now work on older CPUs, too.

Code: Select all

; Experiments In Ray Tracing 2 
; some pretty spheres in realtime (around 70 - 80 FPS on 4GHz Skylake) 
; PB 5.42, done on Windows7 64bit
; Hades, 2016

CompilerIf #PB_Compiler_Thread = 0
  MessageRequester( "Sorry...", "Please, set the compiler to threadsafe before compiling" )
  End
CompilerEndIf

EnableExplicit
DisableDebugger

#TileSize = 16

#RenderWidth = 640  ; RenderWidth/Height needs to be a multiple of #TileSize, otherwise you'll crash!
#RenderHeight = 384 ; Window size and render size are independent, but choosing one as the multiple of the other gives better results.
#WindowWidth = 640
#WindowHeight = 384

#MaxDist = 1000000
#FpError = 0.001

Structure stVector
  x.f
  y.f
  z.f
EndStructure

Structure stMatrix
  m.f[16]
EndStructure  

Structure stColor
  r.f
  g.f
  b.f
EndStructure

Structure stMaterial
  DiffuseColor.stColor
  DiffuseIntensity.f
  SpecularIntensity.f
  Glossiness.f
  AmbientIntensity.f
  Reflectance.f
  Transmittance.f
  IndexOfRefraction.f
  Shader.i
EndStructure
Global NewList Material.stMaterial()

Structure stLight   ; directional light, far away, like the sun
  Direction.stVector
EndStructure
Global Light.stLight

Structure stRay
  Origin.stVector
  Direction.stVector
  length.f
  PrimitiveHit.i
  Inside.i
  Weight.f
EndStructure

Structure stSphere
  Position.stVector
  radius.f
  *Material.stMaterial
EndStructure
Global NewList Sphere.stSphere()
Global Dim SphereX.f(3)
Global Dim SphereY.f(3)
Global Dim SphereZ.f(3)
Global Dim SphereRadius.f(3)
Global Dim *SphereMaterial.stMaterial(3)

Structure stRenderBuffer
  pAddr.i[2]  ; I am using 2 render buffers, so while the last threads are still working on the...
  Width.i     ; ...current frame, the one's that are done can already start with the next one,...
  Height.i
EndStructure
Global RenderBuffer.stRenderBuffer

Structure stCamera
  Pos.stVector
  Matrix.stMatrix
  VScrTL.stVector
  VScrR.stVector
  VScrB.stVector
EndStructure
Global Dim Camera.stCamera(1) ; ...but that means i need two versions for stuff like this, too.

Structure stJob
  x.i
  y.i
  Frame.i
EndStructure
Global NumThreads
Global Dim Job.stJob(0)
Global numJobs
Global JobIdx
Global mutJob
Global mutFrameSync
mutJob = CreateMutex()
mutFrameSync = CreateMutex()

Structure stAABB ; Axis aligned bounding box
  min.f[3]
  max.f[3]
EndStructure

Structure stAABoxHit
  Near.f
  Far.f
EndStructure  

Structure stBVH
  min.f[3]
  max.f[3]
  Content.i
EndStructure

Structure stBVHPrimitive
  AABB.stAABB
  Centre.f[3]
EndStructure  
Global Dim BVHPrimitive.stBVHPrimitive(0)
Global Dim BVHPrimIdx(0)
Global FreeNode.i


Structure stMain
  Running.i
  Frame.i
  DisplList.i
  Material_Default.i
  Shader_Default.i
  Material_Background.i
EndStructure
Global Main.stMain


Macro ColorSet(Color, cr, cg, cb)
  Color\r = cr
  Color\g = cg
  Color\b = cb
EndMacro

Macro VectorSet(Vector, VX, VY, VZ)
  Vector\x = VX
  Vector\y = VY
  Vector\z = VZ
EndMacro

Macro Normalize(Vector, rLength)
  rLength =  1.0 / Sqr(Vector\x * Vector\x + Vector\y * Vector\y + Vector\z * Vector\z)
  Vector\x * rLength
  Vector\y * rLength
  Vector\z * rLength
EndMacro  

Declare Shade(*Ray.stRay, *Color.stColor, Depth.i, ShadowRayFlag.i = #False)

Procedure.i RTMaterialCreate(r.f = 1.0, g.f = 1.0, b.f = 1.0, Shader.i = 0)
  If Shader = 0
    Shader = Main\Shader_Default
  EndIf
  AddElement(Material())
  ColorSet(Material()\DiffuseColor, r, g, b)
  Material()\SpecularIntensity = 0.3
  Material()\Glossiness = 5
  Material()\AmbientIntensity = 0.1
  Material()\Reflectance = 0.0
  Material()\Transmittance = 0.0
  Material()\IndexOfRefraction = 1.0
  Material()\DiffuseIntensity = 1.0 - 1.0 * (Material()\Reflectance + Material()\Transmittance)
  Material()\Shader = Shader
  ProcedureReturn @Material()
EndProcedure

Procedure SphereCreate(x.f, y.f, z.f, radius.f, *Material.stMaterial)
  AddElement(Sphere())
  VectorSet(Sphere()\Position, x, y, z)
  Sphere()\radius = radius
  Sphere()\Material = *Material
EndProcedure

Procedure DisplayBufferInit(BufferWidth.i, BufferHeight.i, VisibleWidth.i, VisibleHeight.i)
  Static Tex.i
  glShadeModel_(#GL_FLAT)
  glEnable_(#GL_TEXTURE_2D)
  glBindTexture_(#GL_TEXTURE_2D, @Tex)
  glTexParameteri_(#GL_TEXTURE_2D, #GL_TEXTURE_MAG_FILTER, #GL_LINEAR)
  glTexParameteri_(#GL_TEXTURE_2D, #GL_TEXTURE_MIN_FILTER, #GL_LINEAR)
  glTexParameteri_(#GL_TEXTURE_2D, #GL_TEXTURE_WRAP_S, #GL_CLAMP)
  glTexParameteri_(#GL_TEXTURE_2D, #GL_TEXTURE_WRAP_T, #GL_CLAMP)
  Main\DisplList = glGenLists_(1)
  glNewList_(Main\DisplList, #GL_COMPILE)
  glBegin_(#GL_QUADS)
  glTexCoord2f_(0.0, 0.0)
  glVertex2i_(0,0)
  glTexCoord2f_(0.0, 1.0)
  glVertex2i_(0, VisibleHeight)
  glTexCoord2f_(1.0, 1.0)
  glVertex2i_(VisibleWidth, VisibleHeight)
  glTexCoord2f_(1.0, 0.0)
  glVertex2i_(VisibleWidth, 0)
  glEnd_()
  glEndList_()
EndProcedure

Procedure glWindowOpen(Width.i, Height.i, Name.s = "Experiments In Ray Tracing 2", flags.i = #PB_Window_SystemMenu | #PB_Window_MinimizeGadget | #PB_Window_ScreenCentered)
  OpenWindow(0, 0, 0, Width, Height, Name, flags)
  OpenGLGadget(0, 0, 0, Width, Height, #PB_OpenGL_NoFlipSynchronization | #PB_OpenGL_NoDepthBuffer)
  SetGadgetAttribute(0, #PB_OpenGL_FlipBuffers, #True)
  glOrtho_(0.0, Width, Height, 0.0, 0.0, 1.0)
EndProcedure

Procedure DisplayBufferShow(pBuffer.i)
  glEnable_(#GL_TEXTURE_2D)
  glTexImage2D_(#GL_TEXTURE_2D, 0, #GL_RGBA, #RenderWidth, #RenderHeight, 0, #GL_BGRA_EXT, #GL_UNSIGNED_BYTE, pBuffer)
  glCallList_(Main\DisplList)
  SetGadgetAttribute(0, #PB_OpenGL_FlipBuffers, #True)
EndProcedure

Procedure Intersect_Spheres(*Ray.stRay)
  Protected MaxD.f, t1.f, t2.f, idx1, idx2, *SphereX.float, *SphereY.float, *SphereZ.float, *SphereRadius.float
  *SphereX.float = @SphereX(0)
  *SphereY.float = @SphereY(0)
  *SphereZ.float = @SphereZ(0)
  *SphereRadius.float = @SphereRadius(0)
  MaxD = #MaxDist
  !mov edx, [p.p_Ray]
  !movss xmm0, [edx]        ; *Ray\Origin\x
  !shufps xmm0, xmm0, 0
  !movss xmm1, [edx + 4]    ; *Ray\Origin\y
  !shufps xmm1, xmm1, 0
  !movss xmm2, [edx + 8]    ; *Ray\Origin\z
  !shufps xmm2, xmm2, 0

  !mov eax, [p.p_SphereX]
  !subps xmm0, [eax]              ; RelPos\x = *Ray\Origin\x - *Sphere\Position\x
  !mov ecx, [p.p_SphereY]
  !subps xmm1, [ecx]              ; RelPos\y = *Ray\Origin\y - *Sphere\Position\y
  !mov eax, [p.p_SphereZ]
  !subps xmm2, [eax]              ; RelPos\z = *Ray\Origin\z - *Sphere\Position\z
  !movss xmm3, [edx + 12]   ; *Ray\Direction\x
  !shufps xmm3, xmm3, 0
  !movss xmm4, [edx + 16]   ; *Ray\Direction\y
  !shufps xmm4, xmm4, 0
  !movss xmm5, [edx + 20]   ; *Ray\Direction\z
  !shufps xmm5, xmm5, 0  
  ;b = *Ray\Direction\x * RelPos\x + *Ray\Direction\y * RelPos\y + *Ray\Direction\z * RelPos\z
  !mulps xmm3, xmm0               ; *Ray\Direction\x * RelPos\x  
  !mulps xmm4, xmm1               ; *Ray\Direction\y * RelPos\y  
  !mulps xmm5, xmm2               ; *Ray\Direction\z * RelPos\z
  
  !addps xmm3, xmm4 
  !addps xmm3, xmm5
  ;c = (RelPos\x * RelPos\x + RelPos\y * RelPos\y + RelPos\z * RelPos\z) - *Sphere\radius * *Sphere\radius
  !mov edx, [p.p_SphereRadius]
  !movups xmm7, [edx]
  !mulps xmm7, xmm7               ; *Sphere\radius * *Sphere\radius
  
  !mulps xmm0, xmm0               ; RelPos\x * RelPos\x  
  !mulps xmm1, xmm1               ; RelPos\y * RelPos\y  
  !mulps xmm2, xmm2               ; RelPos\z * RelPos\z
  
  !addps xmm0, xmm1
  !addps xmm0, xmm2
  !subps xmm0, xmm7
  
  ;d = b * b - c
  !movaps xmm4, xmm3
  !mulps xmm4, xmm4
  !subps xmm4, xmm0
  ; setting everything below 0.0 to 0.0, to not cause trouble with sqrt
  !xorps xmm6, xmm6               ; all 0.0
  !maxps xmm4, xmm6
  
  !sqrtps xmm4, xmm4              ; d = sqr(d)
  
  ; b = -b
  !pcmpeqw xmm7, xmm7             ; all bits to 1
  !pslld xmm7, 31                 ; only sign bits set
  !xorps xmm3, xmm7               ; flip sign bits
  
  ; There are two intersections with a sphere possible: t1 = b - d and t2 = b + d
  !movaps xmm1, xmm3
  !subps xmm3, xmm4               ; t1 = b - d
  !addps xmm1, xmm4               ; t2 = b + d
  ; set t1 to #MaxDist if d <= 0.0 or t1 <= 0.0, so no hit will be detected, to avoid conditonals
  !movss xmm7, [p.v_MaxD]         ; load #MaxDist
  !shufps xmm7, xmm7, 0
  !cmpps xmm4, xmm6, 2            ; create mask d <= 0
  !movaps xmm0, xmm3
  !cmpps xmm0, xmm6, 2            ; create mask t1 <= 0
  !orps xmm0, xmm4                ; combine
  
  ;!blendvps xmm3, xmm7, xmm0      ; t = #MaxDist for d <= 0 or t1 <= 0
  !movaps xmm2, xmm0
  !movaps xmm5, xmm3
  !andps xmm7, xmm2
  !andnps xmm2, xmm5
  !orps xmm7, xmm2
  !movaps xmm3, xmm7

  ; find closest hit
  !movaps xmm5, xmm3
  !movaps xmm0, xmm3
  !shufps xmm0, xmm5, 01001110b   ; shuffle to 1032
  !minps xmm5, xmm0
  !movaps xmm2, xmm5
  !shufps xmm5, xmm5, 10110001b   ; shuffle to 2301
  !minps xmm5, xmm2               ; now smallest t1 is in all positions
  !movss [p.v_t1], xmm5           ; store closest hit
  !cmpeqps xmm5, xmm3             ; set mask for smallest t1 in original
  !movmskps eax, xmm5             ; move sign to rax
  !bsf eax, eax                   ; turn sign into index
  !mov [p.v_idx1], eax            ; store index of closest hit
  
  ; same for t2 now:
  !movaps xmm0, xmm1
  !cmpps xmm0, xmm6, 2
  !orps xmm0, xmm4
  ;!blendvps xmm1, xmm7, xmm0
  !movaps xmm2, xmm0
  !movaps xmm5, xmm1
  !andps xmm7, xmm2
  !andnps xmm2, xmm5
  !orps xmm7, xmm2
  !movaps xmm1, xmm7
  
  !movaps xmm5, xmm1
  !movaps xmm0, xmm1
  !shufps xmm0, xmm5, 01001110b
  !minps xmm5, xmm0
  !movaps xmm2, xmm5
  !shufps xmm5, xmm5, 10110001b
  !minps xmm5, xmm2
  !movss [p.v_t2], xmm5
  !cmpeqps xmm5, xmm1
  !movmskps eax, xmm5
  !bsf eax, eax
  !mov [p.v_idx2], eax
  If t1 < *Ray\length
    *Ray\length = t1
    *Ray\PrimitiveHit = idx1
  EndIf  
  If t2 < *Ray\length
    *Ray\length = t2
    *Ray\PrimitiveHit = idx2
    *Ray\inside = 1
  EndIf  
EndProcedure

Procedure Intersect_Spheres1(*Ray.stRay) ; only here for reference
  Protected RelPos.stVector, i, b.f, c.f, d.f, t1.f, t2.f
  FirstElement(Sphere())
  For i = 0 To ListSize(Sphere()) - 1
    VectorSet(RelPos, *Ray\Origin\x - SphereX(i), *Ray\Origin\y - SphereY(i), *Ray\Origin\z - SphereZ(i))
    b = *Ray\Direction\x * RelPos\x + *Ray\Direction\y * RelPos\y + *Ray\Direction\z * RelPos\z
    c = (RelPos\x * RelPos\x + RelPos\y * RelPos\y + RelPos\z * RelPos\z) - SphereRadius(i) * SphereRadius(i)
    d = b * b - c
    If d > 0.0
      d = Sqr(d)
      t1 = (-b - d)
      If t1 > 0.0 
        If t1 < *Ray\length 
          *Ray\length = t1                  
          *Ray\PrimitiveHit = i
        EndIf
      Else  
        t2 = (-b + d)
        If t2 > 0.0 And t2 < *Ray\length ; a hit from inside the sphere
          *Ray\length = t2
          *Ray\PrimitiveHit = i
          *Ray\inside = 1
        EndIf  
      EndIf
    EndIf
  Next
EndProcedure

Procedure.i Shadow(*HitPosition.stVector, *Normal.stVector, *LightDirection.stVector, LightDistance.f)
  Protected ShadowRay.stRay, Col.stColor
  VectorSet(ShadowRay\Origin, *HitPosition\x, *HitPosition\y, *HitPosition\z)
  VectorSet(ShadowRay\Direction, *LightDirection\x, *LightDirection\y, *LightDirection\z)
  ShadowRay\length = LightDistance
  ShadowRay\PrimitiveHit = -1
  
  Intersect_Spheres(ShadowRay)
  Shade(ShadowRay, Col, 1, #True)
  If ShadowRay\PrimitiveHit >= 0
    ProcedureReturn #True
  EndIf
  ProcedureReturn #False
EndProcedure

Procedure GetSphereNormal(SphIdx.i, *Ray.stRay, *HitPosition.stVector, *Normal.stVector)
  Protected rLength.f
  VectorSet(*Normal, *HitPosition\x - SphereX(SphIdx), *HitPosition\y - SphereY(SphIdx), *HitPosition\z - SphereZ(SphIdx))
  Normalize(*Normal, rLength.f)
  If *Normal\x * *Ray\Direction\x + *Normal\y * *Ray\Direction\y + *Normal\z * *Ray\Direction\z > 0.0
    VectorSet(*Normal, -*Normal\x, -*Normal\y, -*Normal\z)
  EndIf
EndProcedure

Procedure GetLightSpecularColor(*Ray.stRay, *Light.stLight, *Normal.stVector, *LightDirection.stVector, LightNormalAngle.f, *LightSpecularColor.stColor)
  Protected Angle.f, ReflectedLightVector.stVector, Intensity.f, *Material.stMaterial
  *Material = *SphereMaterial(*Ray\PrimitiveHit)
  Angle = 2.0 * LightNormalAngle
  VectorSet(ReflectedLightVector, *LightDirection\x - Angle * *Normal\x, *LightDirection\y - Angle * *Normal\y, *LightDirection\z - Angle * *Normal\z)
  Angle = *Ray\Direction\x * ReflectedLightVector\x + *Ray\Direction\y * ReflectedLightVector\y + *Ray\Direction\z * ReflectedLightVector\z
  If Angle > 0.0
    Intensity = *Material\SpecularIntensity * Pow(Angle, *Material\Glossiness)
    ColorSet(*LightSpecularColor, Intensity, Intensity, Intensity)
  EndIf
EndProcedure

Procedure GetLight(*Ray.stRay, *HitPosition.stVector, *Normal.stVector, *LightDiffuseColor.stColor, *LightSpecularColor.stColor)
  Protected LightNormalAngle.f, LightDistance.f, ShadowRayOrigin.stVector, *Material.stMaterial
  *Material = *SphereMaterial(*Ray\PrimitiveHit)
  ColorSet(*LightDiffuseColor, 0.0, 0.0, 0.0)
  ColorSet(*LightSpecularColor, 0.0, 0.0, 0.0)
  LightDistance = #MaxDist
  LightNormalAngle = *Normal\x * Light\Direction\x + *Normal\y * Light\Direction\y + *Normal\z * Light\Direction\z
  ColorSet(*LightDiffuseColor, *Material\AmbientIntensity, *Material\AmbientIntensity, *Material\AmbientIntensity)
  If LightNormalAngle > *Material\AmbientIntensity
    VectorSet(ShadowRayOrigin, *HitPosition\x + #FpError * *Normal\x, *HitPosition\y + #FpError * *Normal\y, *HitPosition\z + #FpError * *Normal\z)
    If Shadow(ShadowRayOrigin, *Normal, Light\Direction, LightDistance) = #False
      ColorSet(*LightDiffuseColor, LightNormalAngle, LightNormalAngle, LightNormalAngle)
      If *SphereMaterial(*Ray\PrimitiveHit)\SpecularIntensity > 0.0
        GetLightSpecularColor(*Ray, Light, *Normal, Light\Direction, LightNormalAngle, *LightSpecularColor)
      EndIf  
    EndIf
  EndIf
EndProcedure

Procedure Shader_Holes(*Ray.stRay, *Material.stMaterial, *Color.stColor, Depth.i, ShadowRayFlag.i)
  Protected HitPosition.stVector, Normal.stVector, LightDiffuseColor.stColor, LightSpecularColor.stColor, x, y, z, Rel.stVector, d.f, TransmittedRay.stRay
  VectorSet(HitPosition, *Ray\Origin\x + *Ray\length * *Ray\Direction\x, *Ray\Origin\y + *Ray\length * *Ray\Direction\y, *Ray\Origin\z + *Ray\length * *Ray\Direction\z)
  GetSphereNormal(*Ray\PrimitiveHit, *Ray, HitPosition, Normal)
  GetLight(*Ray, HitPosition, Normal, LightDiffuseColor, LightSpecularColor)
  VectorSet(Rel, HitPosition\x * 0.1, HitPosition\y * 0.1, HitPosition\z * 0.1)
  x = Rel\x : y = Rel\y : z = Rel\z
  VectorSet(Rel, Rel\x - x, Rel\y - y, Rel\z - z)
  d = 10 * Sqr(Rel\x * Rel\x + Rel\y * Rel\y + Rel\z * Rel\z)
  If d > 4.8
    If ShadowRayFlag
      ProcedureReturn
    EndIf
    If *Ray\inside
      *Color\r + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\g * LightDiffuseColor\r + LightSpecularColor\r)
      *Color\g + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\r * LightDiffuseColor\g + LightSpecularColor\g)
      *Color\b + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\b * LightDiffuseColor\b + LightSpecularColor\b)
    Else  
      *Color\r + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\r * LightDiffuseColor\r + LightSpecularColor\r)
      *Color\g + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\g * LightDiffuseColor\g + LightSpecularColor\g)
      *Color\b + *Ray\Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\b * LightDiffuseColor\b + LightSpecularColor\b)
    EndIf
  Else
    If ShadowRayFlag
      VectorSet(HitPosition, HitPosition\x - #FpError * Normal\x, HitPosition\y - #FpError * Normal\y, HitPosition\z - #FpError * Normal\z)
      VectorSet(Normal, -Normal\x, -Normal\y, -Normal\z)
      *Ray\PrimitiveHit = -1
      If Shadow(HitPosition, Normal, Light\Direction, #MaxDist)
        *Ray\PrimitiveHit = 0
      EndIf  
      ProcedureReturn
    EndIf  
    VectorSet(TransmittedRay\Origin, HitPosition\x - #FpError * Normal\x, HitPosition\y - #FpError * Normal\y, HitPosition\z - #FpError * Normal\z)
    VectorSet(TransmittedRay\Direction, *Ray\Direction\x, *Ray\Direction\y, *Ray\Direction\z)
    TransmittedRay\length = #MaxDist
    TransmittedRay\PrimitiveHit = -1
    TransmittedRay\Weight = *Ray\Weight
    Intersect_Spheres(TransmittedRay)
    Shade(TransmittedRay, *Color, Depth + 1)
  EndIf  
EndProcedure  

Procedure Shader_Tiles(*Ray.stRay, *Material.stMaterial, *Color.stColor, Depth.i, ShadowRayFlag.i)
  Protected HitPosition.stVector, Normal.stVector, LightDiffuseColor.stColor, LightSpecularColor.stColor
  Protected Tile.i, TileX.i, TileY.i, Brightness.f
  If ShadowRayFlag
    ProcedureReturn
  EndIf  
  VectorSet(HitPosition, *Ray\Origin\x + *Ray\length * *Ray\Direction\x, *Ray\Origin\y + *Ray\length * *Ray\Direction\y, *Ray\Origin\z + *Ray\length * *Ray\Direction\z)
  GetSphereNormal(*Ray\PrimitiveHit, *Ray, HitPosition, Normal)
  GetLight(*Ray, HitPosition, Normal, LightDiffuseColor, LightSpecularColor)
  TileX = HitPosition\x
  TileY = HitPosition\z
  Brightness = 0.6
  TileX & 30
  TileY & 30
  If TileX And TileY
    Brightness + 0.4
    *Color\r + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\r + LightSpecularColor\r)
    *Color\g + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\g + LightSpecularColor\g)
    *Color\b + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\b + LightSpecularColor\b)
  Else
    *Color\r + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\r)
    *Color\g + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\g)
    *Color\b + *Ray\Weight * (*Material\DiffuseIntensity * Brightness * LightDiffuseColor\b)
  EndIf
EndProcedure

Procedure Shader_Sky(*Ray.stRay, *Material.stMaterial, *Color.stColor, Depth.i, ShadowRayFlag.i)
  Protected Angle.f, Brightness.f
  If ShadowRayFlag
    ProcedureReturn
  EndIf  
  Angle = 0.15 + *Ray\Direction\y * 0.85
  If Angle < 0.0
    Angle = 0.0
  EndIf
  Angle = 1.0 - Angle
  Angle = Angle * Angle
  Brightness = 0.2 + 0.7 * Angle * Angle
  *Color\r + *Ray\Weight * Brightness
  *Color\g + *Ray\Weight * Brightness
  *Color\b + *Ray\Weight * 0.9
EndProcedure

Procedure Shader_Default(*Ray.stRay, *Material.stMaterial, *Color.stColor, Depth.i, ShadowRayFlag.i)
  Protected HitPosition.stVector, Normal.stVector, LightDiffuseColor.stColor, LightSpecularColor.stColor
  Protected ReflectedRay.stRay, Angle.f, TransmittedRay.stRay, OutsideIoR.f, RelativeIoR.f, SinT2.f, v.f, Weight.f, RayNormalAngle.f, Fresnel.f
  If ShadowRayFlag
    ProcedureReturn
  EndIf  
  VectorSet(HitPosition, *Ray\Origin\x + *Ray\length * *Ray\Direction\x, *Ray\Origin\y + *Ray\length * *Ray\Direction\y, *Ray\Origin\z + *Ray\length * *Ray\Direction\z)
  GetSphereNormal(*Ray\PrimitiveHit, *Ray, HitPosition, Normal)
  GetLight(*Ray, HitPosition, Normal, LightDiffuseColor, LightSpecularColor)
  If Depth < 10
    If *Material\Reflectance > 0.0
      RayNormalAngle = (*Ray\Direction\x * Normal\x + *Ray\Direction\y * Normal\y + *Ray\Direction\z * Normal\z)
      Fresnel = Pow((1.0 + RayNormalAngle), 5.0)
      ReflectedRay\Weight = *Ray\Weight * (*Material\Reflectance + (1 - *Material\Reflectance) * Fresnel)
      If ReflectedRay\Weight > 0.003
        VectorSet(ReflectedRay\Origin, HitPosition\x + #FpError * Normal\x, HitPosition\y + #FpError * Normal\y, HitPosition\z + #FpError * Normal\z)
        Angle = 2.0 * (*Ray\Direction\x * Normal\x + *Ray\Direction\y * Normal\y + *Ray\Direction\z * Normal\z)
        VectorSet(ReflectedRay\Direction, *Ray\Direction\x - Angle * Normal\x, *Ray\Direction\y - Angle * Normal\y, *Ray\Direction\z - Angle * Normal\z)
        ReflectedRay\length = #MaxDist
        ReflectedRay\PrimitiveHit = -1
        ReflectedRay\Inside = *Ray\Inside
        Intersect_Spheres(ReflectedRay)
        Shade(ReflectedRay, *Color, Depth + 1)
      EndIf
    EndIf
    If *Material\Transmittance > 0.0
      TransmittedRay\Weight = *Ray\Weight * *Material\Transmittance * (1 - Fresnel)
      If TransmittedRay\Weight > 0.003
        ; we assume that outside of the refractive object is air or vacuum (index of refraction = 1.0)
        OutsideIoR = 1.0
        If *Ray\Inside
          RelativeIoR = *Material\IndexOfRefraction / OutsideIoR
          TransmittedRay\Inside = #False  
        Else
          RelativeIoR = OutsideIoR / *Material\IndexOfRefraction
          TransmittedRay\Inside = #True
        EndIf
        ; we have to make sure the origin of the ray is below the surface
        VectorSet(TransmittedRay\Origin, HitPosition\x - #FpError * Normal\x, HitPosition\y - #FpError * Normal\y, HitPosition\z - #FpError * Normal\z)
        Angle = *Ray\Direction\x * Normal\x + *Ray\Direction\y * Normal\y + *Ray\Direction\z * Normal\z
        SinT2 = RelativeIoR * RelativeIoR * (1.0 - Angle * Angle)
        If SinT2 <= 1.0 ; if there is no total internal reflection:
          v = RelativeIoR * Angle + Sqr(1.0 - SinT2)
          VectorSet(TransmittedRay\Direction, RelativeIoR * *Ray\Direction\x - v * Normal\x, RelativeIoR * *Ray\Direction\y - v * Normal\y, RelativeIoR * *Ray\Direction\z - v * Normal\z)
          TransmittedRay\length = #MaxDist
          TransmittedRay\PrimitiveHit = -1
          Intersect_Spheres(TransmittedRay)
          Shade(TransmittedRay, *Color, Depth + 1)
        EndIf
      EndIf
    EndIf
  EndIf
  
  Weight = *Ray\Weight * (1 - (*Material\Reflectance + *Material\Transmittance)) * (1 - Fresnel)
  *Color\r + Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\r * LightDiffuseColor\r + LightSpecularColor\r)
  *Color\g + Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\g * LightDiffuseColor\g + LightSpecularColor\g)
  *Color\b + Weight * (*Material\DiffuseIntensity * *Material\DiffuseColor\b * LightDiffuseColor\b + LightSpecularColor\b)
EndProcedure  

Procedure Shade(*Ray.stRay, *Color.stColor, Depth.i, ShadowRayFlag.i = #False)
  Protected *Material.stMaterial
  If *Ray\PrimitiveHit >= 0
    *Material = *SphereMaterial(*Ray\PrimitiveHit)
  Else
    *Material = Main\Material_Background
  EndIf
  CallFunctionFast(*Material\Shader, *Ray, *Material, *Color, Depth, ShadowRayFlag.i)
EndProcedure

Procedure RenderTile(TileX.i, TileY.i, Frame.i)
  Protected *Buffer.Long, *Material.stMaterial, x.i, y.i, i.i, rLength.f, Ray.stRay, Col.stColor, OffX.stVector, OffY.stVector, RayLS.stVector
  Protected HitPosition.stVector, Normal.stVector, LightNormalAngle.f, Angle.f, ReflectedLightVector.stVector, SpecIntensity.f, Brightness.f
  OffX\x = (Camera(Frame)\VScrR\x - Camera(Frame)\VScrTL\x) / #RenderWidth
  OffX\y = (Camera(Frame)\VScrR\y - Camera(Frame)\VScrTL\y) / #RenderWidth
  OffX\z = (Camera(Frame)\VScrR\z - Camera(Frame)\VScrTL\z) / #RenderWidth
  OffY\x = (Camera(Frame)\VScrB\x - Camera(Frame)\VScrTL\x) / #RenderHeight
  OffY\y = (Camera(Frame)\VScrB\y - Camera(Frame)\VScrTL\y) / #RenderHeight
  OffY\z = (Camera(Frame)\VScrB\z - Camera(Frame)\VScrTL\z) / #RenderHeight
  VectorSet(Ray\Origin, Camera(Frame)\Pos\x, Camera(Frame)\Pos\y, Camera(Frame)\Pos\z)
  For y = TileY To TileY + #TileSize - 1 
    *Buffer = RenderBuffer\pAddr[Frame] + TileX * 4 + y * #RenderWidth * 4
    VectorSet(RayLS, 0.5 + Camera(Frame)\VScrTL\x + y * OffY\x, 0.5 + Camera(Frame)\VScrTL\y + y * OffY\y, Camera(Frame)\VScrTL\z + y * OffY\z)
    For x = TileX To TileX + #TileSize - 1
      VectorSet(Ray\Direction, RayLS\x + x * OffX\x, RayLS\y + x * OffX\y, RayLS\z + x * OffX\z)
      Normalize(Ray\Direction, rLength)
      Ray\PrimitiveHit = -1
      Ray\length = #MaxDist
      Ray\Weight = 1.0
      Ray\Inside = 0
      ColorSet(Col, 0, 0, 0)
      Intersect_Spheres(Ray)
      Shade(Ray, Col, 1)
      *Buffer\l = RGB(Col\b * 255, Col\g * 255, Col\r * 255)
      *Buffer + 4
    Next
  Next
EndProcedure  

Procedure WorkerThread(ThreadIdx)
  Protected iJob, pTravStack.i
  While Main\Running
    iJob = -1
    LockMutex(mutJob)
    If JobIdx < numJobs
      If JobIdx + 1 = numJobs
        LockMutex(mutFrameSync)
      EndIf  
      iJob = JobIdx
      JobIdx + 1    
    EndIf  
    UnlockMutex(mutJob)
    If iJob > -1
      RenderTile(Job(iJob)\x, Job(iJob)\y, Job(iJob)\Frame)
      Job(iJob)\Frame = 1 - Job(iJob)\Frame
      If iJob + 1 = numJobs
        UnlockMutex(mutFrameSync)
      EndIf  
    Else  
      Delay(0)
    EndIf  
  Wend
EndProcedure

Procedure SetUpJobs(RenderWidth.i, RenderHeight.i)
  Protected x, y, i
  ReDim Job((RenderHeight * RenderWidth) / (#TileSize * #TileSize))
  For y = 0 To RenderHeight - 1 Step #TileSize
    For x = 0 To RenderWidth - 1 Step #TileSize
      Job(i)\x = x
      Job(i)\y = y
      i + 1
    Next x  
  Next y  
  numJobs = i
  JobIdx = i
EndProcedure

Procedure CreateThreads()
  Protected i
  NumThreads = CountCPUs(#PB_System_ProcessCPUs)
  For i = 0 To NumThreads - 1
    CreateThread(@WorkerThread(), i)
  Next
EndProcedure

Procedure SceneCreate()
  Protected *Material.stMaterial, i, x, y, z, rLength.f
  *Material = RTMaterialCreate(0.0, 0.0, 0.0, @Shader_Tiles())
  *Material\SpecularIntensity = 0.3
  *Material\Glossiness = 20
  SphereCreate(0, -5025, 0, 5000, *Material)
  
  *Material = RTMaterialCreate(0.7, 0.2, 0.1, @Shader_Holes())
  *Material\SpecularIntensity = 0.2
  *Material\Glossiness = 5
  SphereCreate(-40, -6, 40, 20, *Material)
  
  *Material = RTMaterialCreate(0.0, 0.2, 0.4)
  *Material\Reflectance = 0.01
  *Material\Transmittance = 0.95
  *Material\IndexOfRefraction = 1.5
  *Material\SpecularIntensity = 0.0
  SphereCreate(14, 0, 0, 25, *Material)
  
  *Material = RTMaterialCreate(0.0, 0.2, 0.4)
  *Material\Reflectance = 0.80
  *Material\SpecularIntensity = 0.0
  SphereCreate(5, 5, 50, 15, *Material)
  
  VectorSet(Light\Direction, 1, 1, -1)
  Normalize(Light\Direction, rLength)
  
  FirstElement(Sphere())
  For i = 0 To ListSize(Sphere()) - 1       ; we need to turn the primitive data
    SphereX(i) = Sphere()\Position\x        ; into a SIMD (SSE) friendly format
    SphereY(i) = Sphere()\Position\y        ; Lists and multithreading don't play nice, anyway
    SphereZ(i) = Sphere()\Position\z
    SphereRadius(i) = Sphere()\radius
    *SphereMaterial(i) = Sphere()\Material
    NextElement(Sphere())
  Next i
EndProcedure

Procedure MatrixIdentity(*Matrix.stMatrix)
  *Matrix\m[0] = 1.0 : *Matrix\m[1] = 0.0 : *Matrix\m[2] = 0.0 : *Matrix\m[3] = 0.0
  *Matrix\m[4] = 0.0 : *Matrix\m[5] = 1.0 : *Matrix\m[6] = 0.0 : *Matrix\m[7] = 0.0
  *Matrix\m[8] = 0.0 : *Matrix\m[9] = 0.0 : *Matrix\m[10] = 1.0 : *Matrix\m[11] = 0.0
  *Matrix\m[12] = 0.0 : *Matrix\m[13] = 0.0 : *Matrix\m[14] = 0.0 : *Matrix\m[15] = 1.0
EndProcedure

Procedure MatrixRotateX(*Matrix.stMatrix, Angle.f)
  Protected Help.f, tsin.f, tcos.f
  tsin = Sin(Angle) 
  tcos = Cos(Angle) 
  
  Help = *Matrix\m[4] * tcos + *Matrix\m[8] * tsin
  *Matrix\m[8] = *Matrix\m[4] * -tsin + *Matrix\m[8] * tcos
  *Matrix\m[4] = Help 
  
  Help = *Matrix\m[5] * tcos + *Matrix\m[9] * tsin
  *Matrix\m[9] = *Matrix\m[5] * -tsin + *Matrix\m[9] * tcos
  *Matrix\m[5] = Help 
  
  Help = *Matrix\m[6] * tcos + *Matrix\m[10] * tsin
  *Matrix\m[10] = *Matrix\m[6] * -tsin + *Matrix\m[10] * tcos
  *Matrix\m[6] = Help 
EndProcedure

Procedure MatrixRotateY(*Matrix.stMatrix, Angle.f)
  Protected Help.f, tsin.f, tcos.f
  tsin = Sin(Angle) 
  tcos = Cos(Angle) 
  
  Help = *Matrix\m[0] * tcos - *Matrix\m[8] * tsin
  *Matrix\m[8] = *Matrix\m[0] * tsin + *Matrix\m[8] * tcos
  *Matrix\m[0] = Help 
  
  Help = *Matrix\m[1] * tcos - *Matrix\m[9] * tsin
  *Matrix\m[9] = *Matrix\m[1] * tsin + *Matrix\m[9] * tcos
  *Matrix\m[1] = Help 
  
  Help = *Matrix\m[2] * tcos - *Matrix\m[10] * tsin
  *Matrix\m[10] = *Matrix\m[2] * tsin + *Matrix\m[10] * tcos
  *Matrix\m[2] = Help 
EndProcedure

Procedure MatrixVectorRotation(*out.stVector, *Vector.stVector, *Matrix.stMatrix)
  Protected TempX.f, TempY.f
  TempX = *Matrix\m[0] * *Vector\x + *Matrix\m[1] * *Vector\y + *Matrix\m[2] * *Vector\z 
  TempY = *Matrix\m[4] * *Vector\x + *Matrix\m[5] * *Vector\y + *Matrix\m[6] * *Vector\z
  *out\z = *Matrix\m[8] * *Vector\x + *Matrix\m[9] * *Vector\y + *Matrix\m[10] * *Vector\z
  *out\x = TempX
  *out\y = TempY
EndProcedure 

Procedure Animate(Angle.f)
  Protected Help.f, tsin.f, tcos.f
  tsin = Sin(Angle) 
  tcos = Cos(Angle) 
  Camera(Main\Frame)\Pos\x = -150 * tsin
  Camera(Main\Frame)\Pos\y = 30
  Camera(Main\Frame)\Pos\z = -150 * tcos
  MatrixIdentity(Camera(Main\Frame)\Matrix)
  MatrixRotateX(Camera(Main\Frame)\Matrix, -0.3)
  MatrixRotateY(Camera(Main\Frame)\Matrix, -Angle)
  VectorSet(Camera(Main\Frame)\VScrTL, -0.5 * #RenderWidth, 0.5 * #RenderHeight, #RenderWidth)          ; instead of rotating every single ray,
  MatrixVectorRotation(Camera(Main\Frame)\VScrTL, Camera(Main\Frame)\VScrTL, Camera(Main\Frame)\Matrix) ; I'm using the camera matrix to rotate
  VectorSet(Camera(Main\Frame)\VScrB, -0.5 * #RenderWidth, -0.5 * #RenderHeight, #RenderWidth)          ; a virtual screen, and interpolate the
  MatrixVectorRotation(Camera(Main\Frame)\VScrB, Camera(Main\Frame)\VScrB, Camera(Main\Frame)\Matrix)   ; rays over it, when creating them.
  VectorSet(Camera(Main\Frame)\VScrR, 0.5 * #RenderWidth, 0.5 * #RenderHeight, #RenderWidth)            ; (just to gain some speed)
  MatrixVectorRotation(Camera(Main\Frame)\VScrR, Camera(Main\Frame)\VScrR, Camera(Main\Frame)\Matrix)
  SphereY(3) = -11 + Abs(15 * Sin(Angle * 5)) ; Bounce
EndProcedure

Procedure Init()
  glWindowOpen(#WindowWidth, #WindowHeight, "Experiments In Raytracing 2")
  RenderBuffer\pAddr[0] = AllocateMemory(#RenderWidth * #RenderHeight * 4)
  RenderBuffer\pAddr[1] = AllocateMemory(#RenderWidth * #RenderHeight * 4)
  RenderBuffer\Width = #RenderWidth
  RenderBuffer\Height = #RenderHeight
  DisplayBufferInit(#RenderWidth, #RenderHeight, #WindowWidth, #WindowHeight)
  Main\Shader_Default = @Shader_Default()
  Main\Material_Default = RTMaterialCreate(1.0, 1.0, 1.0)
  Main\Material_Background = RTMaterialCreate(0.0, 0.0, 0.0, @Shader_Sky())
  SceneCreate()
  SetUpJobs(#RenderWidth, #RenderHeight)
  Main\Running = #True
  CreateThreads()
EndProcedure

Procedure Main()
  Protected time.i, AvTime.f, NewTime.f, Event.i, Frame.i, Angle.f, Speed.f
  Init()
  While Main\Running
    If JobIdx >= numJobs
      LockMutex(mutFrameSync)
      UnlockMutex(mutFrameSync)
      DisplayBufferShow(RenderBuffer\pAddr[Main\Frame])
      JobIdx = 0
      Main\Frame = 1 - Main\Frame
      NewTime = ElapsedMilliseconds() - time
      time = ElapsedMilliseconds()
      If Abs(NewTime - AvTime) > NewTime * 0.3
        AvTime = NewTime
      Else  
        AvTime = AvTime * 0.95 + NewTime * 0.05
      EndIf  
      speed = 0.0003 * AvTime
      If speed > 0.01
        speed = 0.01
      EndIf  
      Animate(Angle) : Angle + speed
      SetWindowTitle(0, "Experiments In Raytracing 2 / FPS: " + Str(1000 / AvTime))
    EndIf  
    Delay(0)
    Event = WindowEvent()
    If Event = #PB_Event_CloseWindow
      Main\Running = #False
    EndIf  
  Wend 
  Delay(10)
EndProcedure

Main()

[DK_PETER]

Very cool.

[PureLust]

Wow ... very nice.

(Seems to be 64bit-only. Creates an IMA at Line 235 if I run it in x86-mode.)

[Hades]

Thanks!

And @PureLust: thanks for the report, seem to have been an alignment problem. I've posted a new version that works here with 64bit and x86 compiler.

[applePi]

Thanks Hades for this gift, even it is too much advanced for me but certainly some users will find it informative today or tomorrow .
it runs instantaneous on my system winxp/32. it is good you have used the OpenGLGadget since it simplify the initiating and preparing the environment for the opengl programs. i see there is (pfd.PIXELFORMATDESCRIPTOR, pixformat.i) from the pre OpenGLGadget era (i think). is it useful in the context of this demo since i have deleted it and the demo runs okay as before as i can see.

[Hades]

I am still very much in the process of (re-)learning everything, and a lot of what I am doing right now feels more like scavenging and recycling of my old code, than actually coding, so there is a lot of stuff that I missing, doing in a awkward manner, or simply wrong. Getting feedback helps, so, many thanks!

[djes]

Terrific

[bbanelli]

Code: Select all

; some pretty spheres in realtime (around 70 - 80 FPS on 4GHz Skylake)

Is it possible that 3.1 GHz Haswell (4440) gets only 33-35 FPS?

[Hades]

@bbanelli

That's very strange. On the Skylake worst case is 66-67, best case is around 85, but most of the time it is between 70 and 80. I, too, would have expected that a 3.1 GHz Haswell would be much closer.
Does it run on 8 threads for you? My CPU usage is 96-100%, is that similar for you, too?

Edit: Another thought: How is your CPU usage without running the program? Is there maybe something running in the background, that eats a lot of CPU-cycles?

[Bananenfreak]

Heyho,

I would like to test your code, but it the Compiler throws out an error on line 1 "; Experiments In Ray Tracing 2 ":

Illegal instruction. (Excecution of binary data?)

I´m using the german Version, so this error message is translated by myself. The original message is:

[ERROR] Illegale Anweisung. (Ausführen von binären Daten?)

I'm using PB 5.42 x64, Win 7 Home Premium. Do you have an idea what could be the reason for this error?

[PureLust]

I get 19-24 FPS on my i7-4500U @ 1,8GHz (max 2.4GHz) - all 4 logical CPUs are used.
Still looks great.

I'm using PB 4.2 x64, Win 7 Home Premium. Do you have an idea what could be the reason for this error?

Maybe an Update to PB 5.4 will help?

[bbanelli]

@bbanelli

That's very strange. On the Skylake worst case is 66-67, best case is around 85, but most of the time it is between 70 and 80. I, too, would have expected that a 3.1 GHz Haswell would be much closer.
Does it run on 8 threads for you? My CPU usage is 96-100%, is that similar for you, too?

Edit: Another thought: How is your CPU usage without running the program? Is there maybe something running in the background, that eats a lot of CPU-cycles?

Everything is OK in idle.

http://imgur.com/v76g2LP

Perhaps HT does magic on your Skylake? Anyone with i7 Haswell to confirm?

[Bananenfreak]

I'm using PB 4.2 x64, Win 7 Home Premium. Do you have an idea what could be the reason for this error?

Maybe an Update to PB 5.4 will help?

Upps I meant 5.42

[Hades]

@bbanelli

I've tried to calculate how much FPS you should get on several ways, and every time 75 FPS for me translate to 48-50 for you. When I do the same for PureLust I get 21 FPS, so he is even slightly above expectations. I'm at a loss. I assume your computer is doing fine in benchmarks like cinebench?

@Bananenfreak

What's your CPU? You could try to comment out the Intersect_Spheres() procedure, and remove the '1' from Intersect_Spheres1(). I will try to find the problem from here.

@all

Many thanks for the reports and feedback, much appreciated.

Edit: added a link to 4 .exe's (SSE/PB/32/64Bit) in the first post, in case your compiler is the problem.

[Bananenfreak]

Hi Hades,

my CPU is a AMD Phenom II X4 B55 (Unlocked CPU with formerly 2 cores) @3.6GHz. Equal to X4 955 or X4 965.
Ah, I found it. My processor supports SSE 4a, but not 4.1.
Ok, the other procedure makes it. 7 FPS