31:T480,md` Grid acceleration structure is simplest of spatial subdivision-based acceleration technique. As in many other in its category, the main source of acceleration in this technique is derived by restricting the number of objects the ray intersection once we find that the ray intersects the bounding volume of the scene. The method assumes the bounding volume to be a AABB. It uniformly subdivides the bounding box into a uniform 3D grid, and divides the object by the grid element they spatially belong to. If an object overlaps multiple grid cells, then the object is duplicated in the subdivisions corresponding to those grid cells. This process of grid construction is normally done in a preprocessing step. At the ray tracing time, the grid cells along the direction of the ray are first located in front to back order. The ray is intersected with the group of objects belonging to the grid cell as the cell appears along the ray. The process is continued till an intersection is found, or the ray exits the grid cell. Because of the uniform grid structure, the acceleration technique is also known as uniform spatial subdivision (USS) technique.`32:Tff8,md` Grid construction is a two step process: find the grid resolution, and partition the objects into grid cells (duplicate as necessary). A simple way to choose the resolution will be to decide on the average number of objects per grid cell. This number tells us the average number of ray-intersection tests that must be carried out when the ray intersects a grid cell. So we would like to have this number as small as possible, say n. If we try to make the grid cell to be as close as possible to a cube, then we have a simple algorithm to compute the number of subdivision along each cube dimentions.

let boxVolume = (xSize × ySize × zSize).
let cellVolume = a³ where a the size of the cubic cell is unknown.
numVoxels = boxVolume/cellVolume = N/n where N is the total objects and n average # of objects for cell
So a = (n×boxVolume/N)^1/3 = n^1/3×(boxVolume/N)^1/3 
xResolution = clamp(round(xSize/a),1,MaxResolution)
yResolution = clamp(round(ySize/a),1,MaxResolution)
zResolution = clamp(round(ySize/a),1,MaxResolution)

The algorithm clamps the number of subdivisions along each axis to a user-specified maximum, and this limits the grid resolution. Choice of resolution affects the overhead associated with the acceleration structure. The higher the grid resolution, the smaller the size of the object partition, hence faster is the ray-object-partition intersection, but the higher is the memory overhead (to hold the cell grid) and larger the number of cells to traverse along the ray. The latter can counterbalance the speed-up derived from the acceleration structure. [PBRT](https://www.pbrt.org/) textbook suggests using 64 as the maximum resolution along each axis. Instead of the box volume, it uses a cube of the maximum box dimension and uses n=3³. Under these assumptions the equation for a simplifies to:

a =3×max(xSize, ySize, zSize)/N^1/3.

Once the resolution is decided, partitioning the bounding box to grid cells is straight forward.
Next, we partition the objects into grid cells. As we mentioned earlier, an object may overlap multiple cells. So instead of trying to find out the exact overlap of the object with the grid cell, we find the grid cell indices of the min and max bounds and add the object to the partition of all those cells in between. This process may lead to unnecessary duplication of the objects, but gives us a simple and robust algorithm.

(i_min,j_min,k_min) = gridIndex(minBound)
(i_max,j_max,k_max) = gridIndex(maxBound)
Iteratively add object to the partitions of all the grid cells (i,j,k)
   where i ∈[i_min..i_max], j ∈[j_min..j_max] and k ∈[k_min..k_max].

Grid index for a point is computed as

gridIndex(P)=((x_p-xMin_AABB)/xSize_AABB,(y_p-yMin_AABB)/ySize_AABB,(z_P-zMin_AABB)/zSize_AABB)

As we know, valid index range in each dimension are ([0..xResolution-1], [0..yResolution-1],[0..zResolution-1]). For the rays entering the bounding box from the maximum end, the above equation will produce an index xResolution/yResolution/zResolution. Even if the ray is entering the volume from the minimum end, because of floating-point precision issue, we may get an index as -1. So we must clamp the indices less than 0 to 0 and indices above the max index to max index. i.e.

gridIndex(P)=(clamp((x-xMin_AABB/xSize_AABB), 0, xResolution-1),
	clamp((y-yMin_AABB)/ySize_AABB, 0, yResolution-1),
	clamp((z-zMin_AABB)/zSize_AABB, 0, zResolution-1))

Because of the possibility of duplication, instead of storing object itself, it may be preferable to store object indices in the partitions. `33:T1267,md` In this step we must locate the grid cells in the order in which they are pierced and then intersect the ray with all the objects in the cell partitions in that order, and continue till we find a nearest intersection or till the ray exits the volume, whichever comes first. Nothing is novel about the intersection step. So we will discuss the ordered cell location step. That means, we must find the index of the first pierced cell and then find the next pieced cell and so on. This step must be cariied for every ray, and hence is an overhead to ray-tracing operation. So we must do it as efficiently as possible. At the same time the algorithm must be robust.

To minimize the grid traversal cost, most grid traversal algorithms use incremental operations to find the next cell after a current cell. So we will use such a method as well. If we know the current cell that the ray has pieced, then the next cell must be adjacent along one of the box axes and we need to increment (or decrement) the appropriate index component. Furthermore, as the grid is divided uniformly along each dimension, if we know the intersection parameter t of the ray with the current cell, then the parameter of intersection with the next cell will be simply a dimension specific constant increment. We must be able to detect if the ray exits the volume. As the ray steps through the cells, by incrementing the cell index, ray exiting the volume will be indicated by the index outside the valid index range in each dimension. So we will set up these constant for each ray before we start traversing the cells along the ray.

Δi = (x_rayDirection≥0)?1:-1
Δx=(x_max-x_min)/xResolution
ΔT_x = Δx/abs(x_rayDirection)
OutofRangeIndex_x = (x_rayDirection≥0)?xResolution:-1;
Δj = (y_rayDirection≥0)?1:-1
Δy=(y_max-y_min)/yResolution
ΔT_y = Δy/abs(y_rayDirection)
OutofRangeIndex_y = (y_rayDirection≥0)?yResolution:-1;
Δk = (z_rayDirection≥0)?1:-1
Δz=(z_max-z_min)/zResolution
ΔT_z = Δz/abs(z_rayDirection)
OutofRangeIndex_z = (z_rayDirection≥0)?zResolution:-1;

Next we will determine the index of the first cell pieced by the ray. We know the minimum and maximum parameters of intersection of the ray with the bounding box. So we can compute the point of minimum intersection. Note that we must take into account the possibility that the ray origin is inside the box.

P_intersect = rayPoint(max(0,t_min))
(i,j,k) = gridIndex(P_intersect)

Now we are ready to carry out iterative cell traversal:

while (!intersect ray with objects inside the cell partition){
	find nearest next cell along the ray.
	if the cell is outside volume bound break.
}

Having set up OutofRangeIndex for each dimension, it is easy to find out whether the cell is outside the bound or not.

outside = (i == OutofRangeIndex_x || j == OutofRangeIndex_y || k == OutofRangeIndex_z)

The next cell along the ray must have the index that differs from the current index by at most one increment in either of the dimensions. The nearest next cell is the one that has the smallest t of intersection with the ray. So we will maintain a t for the current cell and maintain next t along the ray direction for each dimension.

t_current=min(0,t_min)
t_x,next=(x_min+(i+Δi)*Δx-x_O)/abs(x_d)
t_y,next=(y_min+(i+Δj)*Δy-y_O)/abs(y_d)
t_z,next=(z_min+(i+Δk)*Δz-z_O)/abs(z_d)

To find the next cell along the ray we find the axis with minimum next t value and increment the index for that axis.

if (t_x,next≤t_y,next&&t_x,next≤t_y,next){i+=Δi; t_current=t_x,next; t_x,next+=Δt_x}
else if (t_y,next≤t_z,next&&t_y,next≤t_x,next){j+=Δj; t_current=t_y,next; t_y,next+=Δt_y}
else {k+=Δk; t_current=t_z,next; t_z,next+=Δt_z}

`34:Tbf3,{ const svgNode = html` `; let svg = d3.select(svgNode); const graph = svg .append("g") .attr("transform", "translate(" + margin.left + "," + margin.top + ")"); let otherDrawings = graph.append("g"); let startPoint = [0, 0]; let endPoint = [0, 0]; let rayDrawing = graph .append("g") .attr("id", "drawn") .style("stroke-width", 1) .style("stroke", "teal") .append("line") .attr("x1", startPoint[0]) .attr("y1", startPoint[1]) .attr("x2", endPoint[0]) .attr("y2", endPoint[1]) .style("marker-start", "url(#dot)") .style("marker-end", "url(#arrow)") .style("opacity", 0); let pickArea = svg .append("rect") .attr("width", width) .attr("height", height) .style("fill", "black") .style("opacity", 0.1) .style("stroke", "black") .on("mousedown", processMouseDown) .on("mouseup", processMouseUp); yGrid(graph); xGrid(graph); graph.append("g").call(d3.axisTop(xScale)); graph.append("g").call(d3.axisLeft(yScale)); graph .append("g") .attr("transform", `translate(${W},0)`) .call(d3.axisRight(yScale)); graph .append("g") .attr("transform", `translate(0 ,${H})`) .call(d3.axisBottom(xScale)); function processMouseDown(event) { event.preventDefault(); // disable text dragging const v = d3.pointer(event); startPoint[0] = v[0] - margin.left; startPoint[1] = v[1] - margin.top; } function processMouseUp(event) { event.preventDefault(); // disable text dragging const v = d3.pointer(event); //mutable endPoint = d3.pointer(event); endPoint[0] = v[0] - margin.left; endPoint[1] = v[1] - margin.top; if (startPoint[0] == endPoint[0] && startPoint[1] == endPoint[1]) return; rayDrawing .attr("x1", startPoint[0]) .attr("y1", startPoint[1]) .attr("x2", endPoint[0]) .attr("y2", endPoint[1]) .style("opacity", 1); let ray = {}; ray.O = { x: xScale.invert(startPoint[0]), y: yScale.invert(startPoint[1]) }; ray.d = { x: xScale.invert(endPoint[0]) - ray.O.x, y: yScale.invert(endPoint[1]) - ray.O.y }; let intersection = rayIntersectsBox(ray, otherDrawings); if (intersection) { renderEndPoints(otherDrawings,intersection.endPoints); const cells = rayTraverseUSS(ray, intersection); renderCells(otherDrawings,ray, cells); pickArea.raise(); } else{ renderCells(otherDrawings,ray, []); renderEndPoints(otherDrawings,[]); } } return svgNode; }35:T7c9,function rayTraverseUSS(ray, intersection) { //Prepare for X let x = ray.O.x + intersection.tMin * ray.d.x; // Point of ray entry let deltaX = (box.max[0] - box.min[0]) / xSubDivisions; let xSlab = Math.floor((x - box.min[0]) / deltaX); if (xSlab >= xSubDivisions) xSlab = xSubDivisions - 1; else if (xSlab < 0) xSlab = 0; let deltaXSlab = ray.d.x < 0 ? -1 : +1; let tx = intersection.tMin; let deltaTx = deltaX / Math.abs(ray.d.x); let xNext = box.min[0] + (xSlab + (ray.d.x < 0 ? 0 : 1)) * deltaX; let tx_next = (xNext - ray.O.x) / ray.d.x; let xSlabLimit = ray.d.x < 0 ? -1 : xSubDivisions; //Prepare for Y let y = ray.O.y + intersection.tMin * ray.d.y; // Point of ray entry let deltaY = (box.max[1] - box.min[1]) / ySubDivisions; let ySlab = Math.floor((y - box.min[1]) / deltaY); if (ySlab >= ySubDivisions) ySlab = ySubDivisions - 1; else if (ySlab < 0) ySlab = 0; let deltaYSlab = ray.d.y < 0 ? -1 : +1; let ty = intersection.tMin; let deltaTy = deltaY / Math.abs(ray.d.y); let yNext = box.min[1] + (ySlab + (ray.d.y < 0 ? 0 : 1)) * deltaY; let ty_next = (yNext - ray.O.y) / ray.d.y; let ySlabLimit = ray.d.y < 0 ? -1 : ySubDivisions; if (gridType == "xSlab") ty_next = Infinity; else if (gridType == "ySlab") tx_next = Infinity; let cells = []; let t = intersection.tMin; while (true) { cells.push([ xSlab, ySlab, t, Math.min(intersection.tMax, Math.min(tx_next, ty_next)) ]); //intersect primitives from cell; //if (intersection in cell) break; if (tx_next < ty_next) { t = tx_next; if (t > intersection.tMax) break; tx_next = tx_next + deltaTx; xSlab += deltaXSlab; if (xSlab == xSlabLimit) break; } else { t = ty_next; if (t > intersection.tMax) break; ty_next = ty_next + deltaTy; ySlab += deltaYSlab; if (ySlab == ySlabLimit) break; } } cells.deltaX = deltaX; cells.deltaY = deltaY; return cells; }36:T476,renderCells = function(svg,ray,cells) { // render cells with opacity 0 svg .selectAll(".segment") .data(cells) .join('line') .attr("class", "segment") .attr("x1", d => xScale(xOnRay(ray, d[2]))) .attr("y1", d => yScale(yOnRay(ray, d[2]))) .attr("x2", d => xScale(xOnRay(ray, d[3]))) .attr("y2", d => yScale(yOnRay(ray, d[3]))) .style("stroke", "red") .style("stroke-width", 1) .style("marker-start", "url(#dot)") .style("marker-end", "url(#arrow)") .style("opacity", 0); svg .selectAll(".cell") .data(cells) .join("rect") .attr("class", "cell") .attr("x", d => gridType == "ySlab" ? 0 : xScale(box.min[0] + d[0] * cells.deltaX) ) .attr("y", d => gridType == "xSlab" ? 0 : yScale(box.min[1] + (d[1] + 1) * cells.deltaY) ) .attr("width", cellWidth) .attr("height", cellHeight) .style("fill", "teal") .style("stroke", "white") .style("opacity", 0); svg .selectAll(".cell,.segment") .transition() .delay(function(d, i) { return i * 300; }) .style("stroke-width", 2) .style("opacity", 0.8); }37:[] 2f:{"title":"Uniform Grid Acceleration","cells":[{"id":346,"value":"md`# Uniform Grid Acceleration`","pinned":false,"hidden":false,"mode":"ojs"},{"id":343,"value":"$31","pinned":false,"hidden":false,"mode":"ojs"},{"id":355,"value":"md`### Grid Construction`","pinned":true,"hidden":false,"mode":"ojs"},{"id":356,"value":"$32","pinned":true,"hidden":false,"mode":"ojs"},{"id":363,"value":"md`### Grid Traversal`","pinned":false,"hidden":false,"mode":"ojs"},{"id":359,"value":"$33","pinned":false,"hidden":false,"mode":"ojs"},{"id":0,"value":"md`# Ray-UniformGrid Traversal demo in 2D`","pinned":false,"hidden":false,"mode":"ojs"},{"id":321,"value":"viewof gridOptions = columns({\n xSubDivisions: slider({\n min: 2,\n max: 16,\n step: 1,\n value: 1,\n description: \"Choose X Subdivisions\"\n }),\n ySubDivisions: slider({\n min: 2,\n max: 16,\n step: 1,\n value: 1,\n description: \"Choose Y Subdivisions\"\n })\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":8,"value":"viewof gridType = radio({\n options:[\"xSlab\",\"ySlab\",\"xyGrid\"],\n value: \"xSlab\"\n})","pinned":false,"hidden":false,"mode":"ojs"},{"id":2,"value":"$34","pinned":false,"hidden":false,"mode":"ojs"},{"id":333,"value":"xSubDivisions = gridOptions.xSubDivisions","pinned":false,"hidden":false,"mode":"ojs"},{"id":334,"value":"ySubDivisions = gridOptions.ySubDivisions ","pinned":false,"hidden":false,"mode":"ojs"},{"id":182,"value":"box = ({ min: [-1, -1], max: [1, 1] })","pinned":false,"hidden":false,"mode":"ojs"},{"id":173,"value":"W = width - margin.left - margin.right","pinned":false,"hidden":false,"mode":"ojs"},{"id":175,"value":"H = height - margin.top - margin.bottom","pinned":false,"hidden":false,"mode":"ojs"},{"id":199,"value":"cellWidth = gridType == \"ySlab\" ? W : W / xSubDivisions","pinned":false,"hidden":false,"mode":"ojs"},{"id":201,"value":"cellHeight = (gridType == \"xSlab\")? H: H / ySubDivisions","pinned":false,"hidden":false,"mode":"ojs"},{"id":167,"value":"xScale = d3\n .scaleLinear()\n .domain([box.min[0], box.max[0]])\n .range([0, W]);","pinned":false,"hidden":false,"mode":"ojs"},{"id":169,"value":"yScale = d3\n .scaleLinear()\n .domain([box.min[1], box.max[1]])\n .range([H, 0]);\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":210,"value":"yAxis = {\n let axis = d3.axisLeft(yScale);\n if (gridType != \"xSlab\")\n axis\n .ticks(xSubDivisions)\n .tickSizeInner(-W)\n .tickSizeOuter(0);\n return axis;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":212,"value":"xAxis = {\n let axis = d3.axisBottom(xScale);\n if (gridType != \"ySlab\")\n axis\n .ticks(ySubDivisions)\n .tickSizeOuter(0)\n .tickSizeInner(-H)\n .tickSizeOuter(0);\n return axis;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":188,"value":"symbolGenerator = d3.symbol().size(20)","pinned":false,"hidden":false,"mode":"ojs"},{"id":190,"value":"crossGenerator = symbolGenerator.type(d3.symbolStar)","pinned":false,"hidden":false,"mode":"ojs"},{"id":276,"value":"d3.range(xScale(box.min[0]), xScale(box.max[0]), cellWidth)","pinned":true,"hidden":false,"mode":"ojs"},{"id":279,"value":"xScale(box.min[0])","pinned":true,"hidden":false,"mode":"ojs"},{"id":281,"value":"xScale(box.max[0])","pinned":true,"hidden":false,"mode":"ojs"},{"id":249,"value":"xGrid = svg =>\n svg.append(\"g\")\n .attr('class', 'grid')\n .style(\"stroke\", \"black\")\n .style(\"stroke-dasharray\", \"5,5\")\n .style(\"opacity\", 0.2)\n .selectAll('line')\n .data(d3.range(xScale(box.min[0]), xScale(box.max[0]), cellWidth))\n .join('line')\n .attr('x1', d => d)\n .attr('x2', d => d)\n .attr('y1', 0)\n .attr('y2', H)","pinned":false,"hidden":false,"mode":"ojs"},{"id":251,"value":"yGrid = svg =>\n svg\n .append(\"g\")\n .attr('class', 'grid')\n .style(\"stroke\", \"black\")\n .style(\"stroke-dasharray\", \"5,5\")\n .style(\"opacity\", 0.2)\n .selectAll('line')\n .data(d3.range(yScale(box.max[1]), yScale(box.min[1]), cellHeight))\n .join('line')\n .attr('x1', 0)\n .attr('x2', W)\n .attr('y1', d => d)\n .attr('y2', d => d)","pinned":false,"hidden":false,"mode":"ojs"},{"id":287,"value":"d3.range(yScale(box.max[1]), yScale(box.min[1]), cellHeight)","pinned":true,"hidden":false,"mode":"ojs"},{"id":291,"value":"yScale(box.min[1])","pinned":true,"hidden":false,"mode":"ojs"},{"id":293,"value":"yScale(box.max[1])","pinned":true,"hidden":false,"mode":"ojs"},{"id":289,"value":"cellHeight","pinned":true,"hidden":false,"mode":"ojs"},{"id":195,"value":"$35","pinned":false,"hidden":false,"mode":"ojs"},{"id":225,"value":"$36","pinned":false,"hidden":false,"mode":"ojs"},{"id":154,"value":"function rayIntersectsBox(ray) {\n let tMinx = (box.min[0] - ray.O.x) / ray.d.x;\n let tMaxx = (box.max[0] - ray.O.x) / ray.d.x;\n if (tMinx > tMaxx) {\n let dummy = tMinx;\n tMinx = tMaxx;\n tMaxx = dummy;\n }\n let tMiny = (box.min[1] - ray.O.y) / ray.d.y;\n let tMaxy = (box.max[1] - ray.O.y) / ray.d.y;\n if (tMiny > tMaxy) {\n let dummy = tMiny;\n tMiny = tMaxy;\n tMaxy = dummy;\n }\n //console.log(tMiny,tMaxy);\n let tMin = tMiny > tMinx ? tMiny : tMinx;\n let tMax = tMaxy < tMaxx ? tMaxy : tMaxx;\n //console.log(tMin,tMax);\n if (tMin < 0) tMin = 0;\n if (tMin >= tMax) return;\n let entryPoint = pointOnRay(ray, tMin);\n let exitPoint = pointOnRay(ray, tMax);\n //console.log(tMin,tMax);\n //console.log(entryPoint, exitPoint);\n return {\n tMin: tMin,\n tMax: Math.min(1, tMax),\n endPoints: [entryPoint,exitPoint]\n };\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":311,"value":"renderEndPoints = (svg, endPoints) => {\n svg\n .selectAll('.symbol')\n .data(endPoints)\n .join('path')\n .attr(\"class\", \"symbol\")\n .attr('transform', function(d) {\n return `translate(${xScale(d.x)} ,${yScale(d.y)})`;\n })\n .attr('d', crossGenerator);\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":156,"value":" function pointOnRay(ray, t) {\n return {\n x: xOnRay(ray, t),\n y: yOnRay(ray, t)\n };\n }\n","pinned":false,"hidden":false,"mode":"ojs"},{"id":158,"value":"function xOnRay(ray, t) {\n return ray.O.x + t * ray.d.x;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":160,"value":"function yOnRay(ray, t) {\n return ray.O.y + t * ray.d.y;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":26,"value":"margin = ({\n top: 100,\n right: 100,\n bottom: 100,\n left: 100\n })","pinned":false,"hidden":false,"mode":"ojs"},{"id":30,"value":"width = 512","pinned":false,"hidden":false,"mode":"ojs"},{"id":32,"value":"height = 512","pinned":false,"hidden":false,"mode":"ojs"},{"id":11,"value":"md`### External Libraries and Imports`","pinned":false,"hidden":false,"mode":"ojs"},{"id":15,"value":"import { radio, slider } from \"@jashkenas/inputs\"","pinned":false,"hidden":false,"mode":"ojs"},{"id":322,"value":"import {columns} from \"@bcardiff/observable-columns\"","pinned":true,"hidden":false,"mode":"ojs"},{"id":13,"value":"d3 = require(\"d3@6\")","pinned":false,"hidden":false,"mode":"ojs"}],"files":"$W37","stdlib":"1","titling":"auto"}