Opticks : Optical Photon Simulation for Particle Physics with NVIDIA OptiX

Opticks Benefits

Outline

• Context
  • Photomultiplier tubes (PMTs)
  • Neutrino detectors
• Current simulation approach, Geant4
• How NVIDIA OptiX can help
• Opticks : optical photon simulation with NVIDIA OptiX
  • geometry migration
  • porting optical physics
• Visualization
• Validation, performance comparisons
• Summary
• Short video

Photomultiplier Tubes (PMTs)

Photomultiplier Tube Operation

Single photoelectron amplified by ~10 dynodes to yield measurable signal

Super-Kamiokande PMTs 1

© Kamioka Observatory, ICRR(Institute for Cosmic Ray Research), The University of Tokyo.

Super-Kamiokande PMTs 3

Super-Kamiokande Expt

• Neutrino interaction produces charged particles travelling faster than c/n ~ 0.75c
• Cherenkov light cone seen as ring by PMTs

© Kamioka Observatory, ICRR(Institute for Cosmic Ray Research), The University of Tokyo.

Super-Kamiokande PMTs 2

© Kamioka Observatory, ICRR(Institute for Cosmic Ray Research), The University of Tokyo.

Super-K Physics

• First observation of neutrino oscillation
  • upward going ~ 0.5 downward going
  • -> non-zero neutrino mass

• Observed supernova 1987a burst neutrinos
  • 11 neutrinos over 10s
  • from nearby galaxy (~160,000 light years)

Jiangmen Underground Neutrino Observatory (JUNO)

Under construction

• neutrino mass heirarchy determination
  • 20 kTon liquid scintillator sphere
  • 53km from 2 nuclear power plants
  • energy resolution: 3% (at 1 Mev)
  • PMTs: 17k 20", 17-35k 3"

• precise oscillation parameter measurement
• study neutrinos from many sources

Dayabay Reactor Neutrino Expt, Far Site

Daya Bay Far Site 2

Antineutrino Detection

Antineutrinos interact with proton of detector --> positron(e+) and neutron (n).

• e+ annihilates, energy deposited in scintillator re-emitted as flash of light
• delayed 2nd flash when gadolinium captures neutron and gamma energy deposited in scintillator

Distinctive double flash of light signature

Daya Bay Far Site 3

Dayabay Physics

Reactor antineutrino measurements

• Most precise oscillation parameter, theta_13
  • Far site baseline ~ 1.9km
  • Near/Far comparison : Antineutrino Disappearance
• Most precise antineutrino spectrum

Daya Bay PMT Wall Photo 1

© 2010 The Regents of the University of California, through the Lawrence Berkeley National Laboratory

Geant4 : Monte Carlo Simulation Toolkit

Geant4 : Monte Carlo Simulation Toolkit Generality

But, there is a problem...

JPMT Before Contact 2

Ray Traced Realistic Image Synthesis and Optical Photon Simulation

OptiX Glass Sample App

http://on-demand.gputechconf.com/siggraph/2013/presentation/SG3106-Building-Ray-Tracing-Applications-OptiX.pdf

• parallel problems, just different focus
• same limit : geometry intersection

Fast ray tracing can revolutionize optical photon simulation

NVIDIA OptiX Ray Tracing Engine:

• state-of-the-art GPU accelerated intersection
• regular releases: improvements, new GPU tuning
• shared C++/CUDA context eases development
• NVIDIA expertise on efficient GPU/multi-GPU usage
  • ~linear scaling with cores across multiple GPUs

NVIDIA OptiX makes massively parallel intersection accessible

Brief History of GPU Optical Photon Simulation Development

2013(Aug-) [liberate geometry from Geant4]

• Develop G4DAE geometry exporter that writes tesselated COLLADA 3D files, including all material and surface properties.

2014 [integrate geometry/event data with Chroma [1]]

• Integrate G4DAE geometries with Chroma
• Develop runtime event data bridge
• GPU Cerenkov/Scintillation photon generation

2015 [replace Chroma with Opticks]

• Develop Opticks with NVIDIA OptiX ray tracing
• Opticks/Geant4 matched for simple geometries
• Speedup factor of 200x with a mobile GPU

2016 [Opticks validation against Geant4]

• Achieved match for single PMT

[1] http://chroma.bitbucket.org Developed by Stan Seibert, University of Pennsylvania.

Migrating Geometry from Geant4 to Opticks

Translate Geant4 CSG tree solids into Opticks boundaries:

• inner material : self
• outer material : parent
• inner surface : outwards going photons (self to parent)
• outer surface : inwards going photons (parent to self)

Properties as function of wavelength interleaved into float4 texture.

Texture lookup simplifies OptiX programs

// material property lookup from OptiX programs
float nmi = (nm - wavelength_domain.x)/wavelength_domain.z + 0.5f ;
float4 material1 = tex2D(boundary_texture, nmi, line + 0.5f );

float refractive_index = material1.x ;
float absorption_length = material1.y ;
float scattering_length = material1.z ;
float reemission_probability = material1.w ;  // >0 for scintillators

Large Geometry/Event Techniques

Geometry analysed to find instances

JUNO: ~90M --> 0.1M triangles

• 19,000 20" PMTs
• 34,000 3" PMTs

OpenGL/OptiX instancing

Multiple Renderers

• global, full instances, bbox instances
• rapid switching controls

Split simulation from visualization

• compute mode
  • pure OptiX buffers ~7X interop
• visualize by separate event loading

Geometry Modelling : Tesselated vs Analytic Photomultiplier Tubes

Analytic : more realistic, faster, less memory, much more effort

For Dayabay PMT:

• partition CSG solids into 12 single primitive parts (instead of 2928 triangles)
• splitting at geometrical intersections avoids implementing general CSG boolean handling
• geometry provided to OptiX in form of ray intersection code

Tesselated Photomultiplier : unrealistic disco ball effect

Analytic PMTs together with triangulated geometry

Starting from convenient tesselated geometry, implement analytic geometry for critical parts in primary optical path, to allow very close Geant4/Opticks match.

Visualizing An Optical Photon Simulation

Compressed Steps, OpenGL Buffer, written by OptiX

• space for 16 steps per photon
  • position, time, polarization, wavelength

OpenGL Geometry Shader: lines -> points/line_strip

• event time is uniform input
• pair validity from step times
• pair selection from step times relative to input time
• interpolation between steps gives position

Interactive: time scrubbing + history/material selection

// geom.glsl
layout (lines) in;
layout (line_strip, max_vertices = 2) out;
void main ()
{
    vec4 p0 = gl_in[0].gl_Position  ;   // p0.w : step begin time
    vec4 p1 = gl_in[1].gl_Position  ;   // p1.w : step end time
    ...

Compare Opticks/Geant4 Simulations with Simple Lights/Geometries

1M Photons -> Water Sphere (S-Polarized)

0.5M Photons -> Dayabay PMT

• record steps of all photons : position, time, wavelength, polarization vector, material/history codes
• 128 bit per step (2*short4) : snorm/uchar compressed using known domains
• water sphere : yields multiple rainbows depending on reflections

Opticks/Geant4 Rainbow Step Sequence Comparison

Statistically consistent photon histories in the two simulations

• BT/BR: boundary transmit/reflect
• TO/SC/SA: torch/scatter/surface absorb

 64-bit uint  Opticks    Geant4    chi2                                      (tag:5,-5)

        8ccd   819160    819654    0.15  [4 ] TO BT BT SA                    (cross droplet)
         8bd   102087    101615    1.09  [3 ] TO BR SA                       (external reflect)
       8cbcd    61869     61890    0.00  [5 ] TO BT BR BT SA                 (bow 1)
      8cbbcd     9618      9577    0.09  [6 ] TO BT BR BR BT SA              (bow 2)
     8cbbbcd     2604      2687    1.30  [7 ] TO BT BR BR BR BT SA           (bow 3)
    8cbbbbcd     1056      1030    0.32  [8 ] TO BT BR BR BR BR BT SA        (bow 4)
       86ccd     1014      1000    0.10  [5 ] TO BT BT SC SA
   8cbbbbbcd      472       516    1.96  [9 ] TO BT BR BR BR BR BR BT SA     (bow 5)
         86d      498       473    0.64  [3 ] TO SC SA
  bbbbbbbbcd      304       294    0.17  [10] TO BT BR BR BR BR BR BR BR BR  (bow 8+ truncated)
  8cbbbbbbcd      272       247    1.20  [10] TO BT BR BR BR BR BR BR BT SA  (bow 6)
  cbbbbbbbcd      183       161    1.41  [10] TO BT BR BR BR BR BR BR BR BT  (bow 7 truncated)

1M 64 bit uint flag sequences indexed using CUDA Thrust, 0.040 s (sparse histogram sort)

Rainbow Spectrum for 1st six bows

Spectra obtained by selecting photons by internal reflection counts. Colors obtained from spectra of each bin using CIEXYZ weighting functions converted into sRGB/D65 colorspace. Exposures by normalizing to bin with maximum luminance (CIE-Y) of each bow. White lines indicate geometric optics prediction of deviation angle ranges of the visible range 380-780nm. 180-360 degrees signifies exit on same side of droplet as incidence.

1M Rainbow S-Polarized, Comparison Opticks/Geant4

Deviation angle(degrees) of 1M parallel monochromatic photons in disc shaped beam incident on water sphere. Numbered bands are visible range expectations of first 11 rainbows. S-Polarized intersection (E field perpendicular to plane of incidence) arranged by directing polarization radially.

PmtInBox after 1

PmtInBox after 2

PMT Opticks/Geant4 step distribution comparison TO BT [SD]

Good agreement reached, after several fixes: geometry, total internal reflection, group velocity

position(xyz), time(t)

polarization(abc), radius(r)

PMT Opticks/Geant4 step distribution comparison : chi2/ndf

Photon Propagation Times Geant4 cf Opticks

• Opticks > 200X Geant4 with only 384 core mobile GPU[1] (multi-GPU workstation up to 20x more cores)
• photon propagation time will become effectively zero

• Interop uses OpenGL buffers allowing visualization, Compute uses OptiX buffers
• Interop/Compute : perfectly identical results, monitored by digest

[1] MacBook Pro (2013), NVIDIA GeForce GT 750M, 2048 MB, 384 cores

Summary

Opticks enables particle physics to benefit from optical photon simulation taking effectively zero time and zero CPU memory, thanks to massive parallelism made accessible by NVIDIA OptiX.

• The more photons the bigger the overall speedup (99% -> 100x)
• Drastic speedup -> better detector understanding -> greater precision
• Large PMT based neutrino experiments, such as JUNO, can benefit the most

https://bitbucket.org/simoncblyth/opticks

YouTube Video

https://www.youtube.com/watch?v=QzH6y0pKXk4

YouKu Video

http://v.youku.com/v_show/id_XMTUxNjM3MzA4OA==.html

JPMT Wide

http://simoncblyth.bitbucket.io/env/graphics/ggeoview/jpmt-wide_half.png

JPMT After Contact

http://simoncblyth.bitbucket.io/env/graphics/ggeoview/jpmt-after-contact_half.png

JPMT Approach

http://simoncblyth.bitbucket.io/env/graphics/ggeoview/jpmt-approach_half.png

JPMT Arrival

http://simoncblyth.bitbucket.io/env/graphics/ggeoview/jpmt-arrival_half.png

JPMT Inside Outside

http://simoncblyth.bitbucket.io/env/graphics/ggeoview/jpmt-inside-outside_half.png

OptiX Experience

Geometry Issues Fixed/Avoidable

• triangulated geometry leaks when shoot millions of photons at cracks (more a problem for testing than usage)
• bounding boxes that touch geometry cause leaks, slightly enlarging the boxes avoids leak (<1 in 3M)

Use cuRAND in OptiX without large stack

Init in separate CUDA launch, persist cuRAND state, subsequently load into OptiX context : avoids large stack performance hit

Compositing OptiX Ray Trace and OpenGL Rasterized

• best of both worlds : rasterization speed, raytrace accuracy
• OptiX : replace uchar4 alpha with Z clip depth
• OpenGL : fragment shader, set gl_FragDepth

out vec4 frag_colour;
uniform sampler2D ColorTex ;
void main ()
{
   frag_colour = texture(ColorTex, texcoord);
   float depth = frag_colour.w ;
   frag_colour.w = 1.0 ;
   gl_FragDepth = depth  ;
   ...