Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

List intersections #6442

Merged
merged 15 commits into from
Nov 3, 2023
Merged

Conversation

dmitrishastin
Copy link
Contributor

@dmitrishastin dmitrishastin commented Oct 22, 2023

Add "ListIntersections" (C++) / "list_intersections" (python) function that returns a dictionary for all ray intersections present in a raycasting scene.

Type

  • Bug fix (non-breaking change which fixes an issue): Fixes #
  • New feature (non-breaking change which adds functionality). Resolves #
  • Breaking change (fix or feature that would cause existing functionality to not work as expected) Resolves #

Motivation and Context

Flexibility in ray inclusion / exclusion.

Checklist:

  • I have run python util/check_style.py --apply to apply Open3D code style
    to my code.
  • This PR changes Open3D behavior or adds new functionality.
    • Both C++ (Doxygen) and Python (Sphinx / Google style) documentation is
      updated accordingly.
    • I have added or updated C++ and / or Python unit tests OR included test
      results
      (e.g. screenshots or numbers) here.
  • I will follow up and update the code if CI fails.
  • For fork PRs, I have selected Allow edits from maintainers.

Description

Accepts a rays tensor. The returned dictionary contains:

ray_ids
    A tensor with ray IDs. The shape is {..}.
    
geometry_ids
    A tensor with the geometry IDs. The shape is {..}.

primitive_ids
    A tensor with the primitive IDs, which corresponds to the triangle
    index. The shape is {..}.

primitive_uvs 
    A tensor with the barycentric coordinates of
    the closest points within the triangles. The shape is {.., 2}.

t_hit
    A tensor with the distance to the hit. The shape is {..}. 

Runs count_intersections first to pre-allocate results arrays then executes the actual function.

Can be tested with the following code:

import open3d as o3d
import numpy as np

d = o3d.data.MonkeyModel()
mesh = o3d.t.io.read_triangle_mesh(d.path)
p_min = np.min(mesh.vertex['positions'].numpy() - 1, axis=0)
p_max = np.max(mesh.vertex['positions'].numpy() - 1, axis=0)
x = np.arange(p_min[0], p_max[0], .1)
y = np.arange(p_min[1], p_max[1], .1)
xv, yv = np.meshgrid(x, y) 
orig = np.vstack([xv.flatten(), yv.flatten(), np.tile(p_min[2], xv.size)]).T
dest = np.copy(orig)
dest[:, 2] = p_max[2] + 1
rays = np.hstack([orig, dest - orig]).astype('float32') 

scene = o3d.t.geometry.RaycastingScene()
mesh_id = scene.add_triangles(mesh)
lx = scene.list_intersections(rays)

Output:

>>> rays.shape
(540, 6)

>>> lx['ray_ids'].numpy().T
array([[318, 318, 319, 319, 320, 320, 321, 321, 322, 322, 323, 323, 345,
        345, 346, 346, 347, 347, 348, 348, 349, 349, 350, 350, 372, 372,
        373, 373, 374, 374, 375, 375, 376, 376, 377, 377, 399, 399, 400,
        400, 401, 401, 402, 402, 403, 403, 404, 404, 426, 426, 427, 427,
        428, 428, 429, 429, 430, 430, 431, 431, 453, 453, 454, 454, 455,
        455, 456, 456, 457, 457, 480, 480, 481, 481, 482, 482, 483, 483,
        484, 484, 485, 485, 505, 505, 506, 506, 507, 507, 508, 508, 509,
        509, 510, 510, 511, 511, 512, 512, 526, 526, 527, 527, 528, 528,
        529, 529, 530, 530, 531, 531, 532, 532, 533, 533, 534, 534, 535,
        535, 536, 536, 537, 537, 538, 538, 539, 539]], dtype=uint32)

>>> lx['geometry_ids'].numpy().T
array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
        0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=uint32)

>>> lx['primitive_ids'].numpy().T
array([[12297,  1938,  4245,  1912, 12109,  2832, 12088,  2787,  4417,
         9773, 12287,  8901,  4433,  2006,  4235,  1946,  4211, 10834,
        12093,  2959, 12295,  1922, 12283,  8938,  4400,  1983, 12085,
         2015, 12077,  2109, 12068,  9961, 12253,  1998, 12225,  1097,
        12279,  8977,  4203,  1969, 12051,  9948,  4194,  9920,  4385,
         9836,  4359,  1096, 12333,  1111, 12042,  1973,  4175,  2065,
        12029,  9927, 12321,  9825, 12356,  1089,  4787, 10512, 12588,
        10488,  4721,  2036,  4696,  9890,  4704,  2606, 12555, 12406,
        12564, 10452,  4694,  2301,  4677,  2277, 12544,  2574,  4675,
         5650,  5823,  5575,  5813, 13514,  4685,  2715,  4690, 10428,
         4691, 10193,  4667, 10312,  4670,  2692, 12543, 13548, 15393,
         6035,  7524, 13900, 15387,  6283, 15728,  6292, 15733, 13451,
        13648, 13476,  5809,  9054, 13710,  9010,  4648,  1718,  4660,
         1759,  4662, 10392,  4645, 10376, 12512,  1703,  4642,  9006]],
      dtype=uint32)

>>>lx['primitive_uvs'].numpy().T
array([[8.50260928e-02, 6.51837885e-01, 1.03081293e-01, 3.91855270e-01,
        8.05648923e-01, 5.80210924e-01, 1.40363634e-01, 4.72955555e-01,
        3.64858925e-01, 2.04821900e-01, 6.26785994e-01, 7.60398880e-02,
        1.04184061e-01, 2.94283777e-01, 2.93564528e-01, 6.82043910e-01,
        1.17763458e-02, 1.93700671e-01, 3.53171467e-03, 3.03192195e-02,
        5.40782988e-01, 1.92260608e-01, 4.34213966e-01, 6.12139523e-01,
        5.72637856e-01, 6.21795833e-01, 1.55210704e-01, 1.16427779e-01,
        9.15720105e-01, 3.63248348e-01, 8.89678299e-01, 5.30207932e-01,
        6.57920063e-01, 8.82471576e-02, 1.01874083e-01, 1.09900184e-01,
        7.86101043e-01, 1.81164563e-01, 5.11571527e-01, 1.69970412e-02,
        1.00328386e-01, 2.37347394e-01, 2.26802994e-02, 1.81675166e-01,
        5.41030526e-01, 7.16622248e-02, 4.13606241e-02, 2.97267050e-01,
        4.64487553e-01, 1.66320521e-02, 5.10603547e-01, 4.46072996e-01,
        6.01617873e-01, 3.57365131e-01, 5.13987005e-01, 2.44175628e-01,
        1.64774716e-01, 2.06305742e-01, 6.90514445e-01, 1.78633600e-01,
        1.05318859e-01, 2.54442692e-01, 4.06445339e-02, 3.29573601e-01,
        5.06643891e-01, 2.37140119e-01, 1.31217435e-01, 2.17609853e-01,
        3.27986866e-01, 4.65133274e-03, 8.08010638e-01, 3.72975886e-01,
        6.18177712e-01, 5.14547825e-01, 4.61467683e-01, 1.66717678e-01,
        6.89924583e-02, 4.14696962e-01, 2.27413654e-01, 7.33788610e-01,
        4.05685335e-01, 7.58898035e-02, 1.38070121e-01, 5.24565458e-01,
        2.51125902e-01, 7.64557898e-01, 4.52516317e-01, 7.91438341e-01,
        2.20954821e-01, 6.51578248e-01, 5.62553287e-01, 4.47666347e-01,
        7.86239356e-02, 1.15206286e-01, 6.36142194e-01, 3.99974853e-01,
        4.65191007e-01, 2.71531105e-01, 1.17488220e-01, 8.56871128e-01,
        1.96953505e-01, 6.84856653e-01, 4.82201636e-01, 2.80162394e-01,
        1.94526419e-01, 5.61448792e-03, 6.73465431e-02, 4.94458348e-01,
        7.45218247e-02, 9.88216773e-02, 1.32805824e-01, 6.08568132e-01,
        7.81525373e-02, 1.48485616e-01, 3.72431427e-01, 1.60334751e-01,
        7.44871572e-02, 4.66253549e-01, 3.14146519e-01, 4.40233052e-02,
        1.48861185e-01, 3.42773944e-01, 3.97519767e-02, 4.67371553e-01,
        4.99460287e-02, 8.69180620e-01],
       [3.92546691e-02, 1.67678565e-01, 3.65194201e-01, 4.47089612e-01,
        1.26379047e-04, 3.68818521e-01, 1.90832779e-01, 1.35176867e-01,
        3.05290520e-01, 5.80557048e-01, 2.97499180e-01, 5.30879274e-02,
        2.21387267e-01, 6.46204174e-01, 3.79589945e-01, 1.86955437e-01,
        8.65287006e-01, 6.06772840e-01, 4.78246175e-02, 6.14985645e-01,
        3.59579146e-01, 2.66485035e-01, 2.73514956e-01, 3.65829051e-01,
        1.82528138e-01, 3.16096485e-01, 4.36164469e-01, 4.94528979e-01,
        2.34494451e-02, 1.03755906e-01, 6.05465584e-02, 4.63592499e-01,
        6.02329001e-02, 2.00819299e-01, 5.75036347e-01, 8.54499578e-01,
        1.67031765e-01, 2.93776039e-02, 6.69649169e-02, 1.01286829e-01,
        8.44411194e-01, 2.02314824e-01, 2.68521532e-02, 7.42257297e-01,
        3.14291596e-01, 1.53447643e-01, 6.90125883e-01, 6.37456998e-02,
        4.60704058e-01, 7.82134652e-01, 4.86005783e-01, 3.26901466e-01,
        2.96064347e-01, 4.09353852e-01, 3.34959358e-01, 6.70883581e-02,
        4.87725213e-02, 6.34807169e-01, 1.99106425e-01, 5.09882152e-01,
        6.14799976e-01, 1.92902938e-01, 2.93571591e-01, 7.06736818e-02,
        4.36887369e-02, 2.95727044e-01, 8.39247465e-01, 9.20358598e-02,
        5.14076471e-01, 1.55397788e-01, 1.48478355e-02, 9.43628401e-02,
        3.18643987e-01, 3.43717992e-01, 1.53597191e-01, 3.62733036e-01,
        1.27671719e-01, 3.86326134e-01, 5.63073039e-01, 7.00732917e-02,
        5.23127794e-01, 7.82335460e-01, 6.78201020e-01, 2.89844066e-01,
        2.95127273e-01, 3.91148739e-02, 4.72639620e-01, 1.01692388e-02,
        1.92888826e-01, 1.89839691e-01, 7.12942183e-02, 6.05209507e-02,
        8.27150643e-01, 3.16473961e-01, 7.92979747e-02, 5.56479931e-01,
        2.22815737e-01, 1.68949321e-01, 3.75169277e-01, 7.28604645e-02,
        7.43002892e-01, 1.27236247e-01, 1.60706714e-01, 6.42337501e-01,
        4.55084860e-01, 5.32845736e-01, 7.88478732e-01, 3.08783919e-01,
        2.64580131e-01, 7.74997652e-01, 6.83005676e-02, 1.53609201e-01,
        4.40520406e-01, 3.76733273e-01, 4.96786356e-01, 7.34959781e-01,
        3.65802705e-01, 5.00856578e-01, 3.25992048e-01, 1.70454092e-03,
        1.85889587e-01, 2.30168775e-01, 4.06972110e-01, 2.14318082e-01,
        5.39930463e-01, 1.18742652e-01]], dtype=float32)

>>> lx['t_hit'].numpy().T
array([[0.8435714 , 0.93022877, 0.84232   , 0.94535124, 0.83926505,
        0.9475919 , 0.8404622 , 0.949586  , 0.8430659 , 0.9405413 ,
        0.8463084 , 0.921236  , 0.83460045, 0.9428818 , 0.829282  ,
        0.958152  , 0.8238311 , 0.94325334, 0.82550997, 0.94131124,
        0.83107746, 0.95320255, 0.839386  , 0.93259716, 0.8331073 ,
        0.9460807 , 0.8237588 , 0.96223986, 0.81628525, 0.9666485 ,
        0.81862414, 0.9659599 , 0.826654  , 0.95780087, 0.8420907 ,
        0.93435174, 0.8350253 , 0.94551253, 0.8188721 , 0.9639249 ,
        0.80940276, 0.9687127 , 0.8120805 , 0.967958  , 0.8238137 ,
        0.9589816 , 0.8503339 , 0.9307414 , 0.8364902 , 0.9418003 ,
        0.8019855 , 0.9645548 , 0.7877976 , 0.97038263, 0.79091233,
        0.9694287 , 0.81348133, 0.9587221 , 0.8602229 , 0.9216753 ,
        0.83709997, 0.9308639 , 0.6462993 , 0.9639895 , 0.6154006 ,
        0.97064835, 0.6210553 , 0.96957433, 0.68068993, 0.95639205,
        0.55123186, 0.8929472 , 0.5219662 , 0.96521115, 0.50918746,
        0.9916076 , 0.51180327, 0.98705995, 0.5309491 , 0.95183134,
        0.5730049 , 0.80399156, 0.59807265, 0.73571694, 0.516519  ,
        0.84292406, 0.48224548, 0.9099531 , 0.4638537 , 0.97253287,
        0.45487374, 0.99615115, 0.4568556 , 0.9969643 , 0.46949258,
        0.9533076 , 0.49270236, 0.8698291 , 0.5398988 , 0.5576007 ,
        0.5450186 , 0.5949125 , 0.56059927, 0.6183859 , 0.56956923,
        0.6249971 , 0.5710189 , 0.70505244, 0.55278033, 0.8113385 ,
        0.50697887, 0.9186829 , 0.4648054 , 0.9423307 , 0.44438577,
        0.964288  , 0.43005514, 0.97802615, 0.42262608, 0.9722668 ,
        0.42434618, 0.9756199 , 0.43452722, 0.9747778 , 0.45160547,
        0.9560417 ]], dtype=float32)

Visualisation:

entry = o3d.geometry.PointCloud(points = o3d.utility.Vector3dVector(orig))
target = o3d.geometry.PointCloud(points = o3d.utility.Vector3dVector(dest))
correspondence = [(i, i) for i in range(rays.shape[0])]
traj = o3d.geometry.LineSet.create_from_point_cloud_correspondences(entry , target , correspondence)
traj.colors = o3d.utility.Vector3dVector(np.tile([1, 0, 0], [rays.shape[0], 1]))
o3d.visualization.draw([mesh, traj])

# intersection points
v = mesh.vertex['positions'].numpy()
t = mesh.triangle['indices'].numpy()
tidx = lx['primitive_ids'].numpy()
uv = lx['primitive_uvs'].numpy()
w = 1 - np.sum(uv, axis=1)
c = \
  v[t[tidx, 1].flatten(), :] * uv[:, 0][:, None] + \
  v[t[tidx, 2].flatten(), :] * uv[:, 1][:, None] + \
  v[t[tidx, 0].flatten(), :] * w[:, None] 
x = o3d.geometry.PointCloud(points = o3d.utility.Vector3dVector(c))

# visualise
o3d.visualization.draw([mesh, traj, x])

image


This change is Reviewable

@update-docs
Copy link

update-docs bot commented Oct 22, 2023

Thanks for submitting this pull request! The maintainers of this repository would appreciate if you could update the CHANGELOG.md based on your changes.

Copy link
Contributor

@benjaminum benjaminum left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks @dmitrishastin ! I think this is a good contribution. Please have a look at the comments. Other things that are missing is an example inside the python docstring and a unit test for the new function. Let me know if you need help with this.

Reviewable status: 0 of 3 files reviewed, 3 unresolved discussions (waiting on @dmitrishastin)


cpp/open3d/t/geometry/RaycastingScene.cpp line 553 at r1 (raw file):

                           const size_t num_rays,
                           const size_t num_intersections,
                           Eigen::VectorXi cumsum,

Pass cumsum as const reference to avoid copy


cpp/open3d/t/geometry/RaycastingScene.cpp line 855 at r1 (raw file):

    Eigen::Map<Eigen::VectorXi> intersections_vector(
            intersections.GetDataPtr<int>(), num_rays);
    Eigen::Map<Eigen::VectorXi> num_intersections(

Isn't the temporary Tensor destroyed here and the pointer invalid?


cpp/pybind/t/geometry/raycasting_scene.cpp line 203 at r1 (raw file):

    
    ray_ids
        A tensor with ray IDs. The shape is {..}.

I think it would be easier to use the result if the cumsum. It would be similar to here http://www.open3d.org/docs/release/python_api/open3d.ml.tf.layers.RadiusSearch.html where each query can have a different number of neighbors. @dmitrishastin what is your opinion here?

Copy link
Contributor Author

@dmitrishastin dmitrishastin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for reviewing promptly @benjaminum ! I addressed the comments to the best of my ability. I have also included the example inside the python docstring as suggested, and checked it by compiling documentation. I am, however, not clear regarding the unit test as I could not find any unit tests for RaycastingScene - your help would be appreciated. Finally, I have expanded my original submitted function to also output primitive_uvs as this made sense.

Reviewable status: 0 of 3 files reviewed, 3 unresolved discussions (waiting on @benjaminum)


cpp/open3d/t/geometry/RaycastingScene.cpp line 553 at r1 (raw file):

Previously, benjaminum (Benjamin Ummenhofer) wrote…

Pass cumsum as const reference to avoid copy

Done.


cpp/open3d/t/geometry/RaycastingScene.cpp line 855 at r1 (raw file):

Previously, benjaminum (Benjamin Ummenhofer) wrote…

Isn't the temporary Tensor destroyed here and the pointer invalid?

Thanks for spotting! I have changed so it hopefully makes more sense.


cpp/pybind/t/geometry/raycasting_scene.cpp line 203 at r1 (raw file):

Previously, benjaminum (Benjamin Ummenhofer) wrote…

I think it would be easier to use the result if the cumsum. It would be similar to here http://www.open3d.org/docs/release/python_api/open3d.ml.tf.layers.RadiusSearch.html where each query can have a different number of neighbors. @dmitrishastin what is your opinion here?

Apologies - I am running into trouble compiling the ML module on Windows. Would you be able to provide an example of RadiusSearch output for me to be more clear on your proposed alternative? Thank you!

Copy link
Contributor

@benjaminum benjaminum left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The unit tests for RaycastingScene are using the python bindings and are here https://github.com/isl-org/Open3D/blob/master/python/test/t/geometry/test_raycasting_scene.py . Adding primitive_uvs makes sense.

Thanks again and see the new comments!

Reviewable status: 0 of 3 files reviewed, 2 unresolved discussions (waiting on @dmitrishastin)


cpp/open3d/t/geometry/RaycastingScene.cpp line 553 at r1 (raw file):

Previously, dmitrishastin (Dmitri Šastin) wrote…

Done.

The & is missing to avoid the copy -> const Eigen::VectorXi& cumsum


cpp/pybind/t/geometry/raycasting_scene.cpp line 203 at r1 (raw file):

Previously, dmitrishastin (Dmitri Šastin) wrote…

Apologies - I am running into trouble compiling the ML module on Windows. Would you be able to provide an example of RadiusSearch output for me to be more clear on your proposed alternative? Thank you!

The idea is to use the cumulative sum for pointing to the start and end of the intersections for each ray similar to this example https://www.tensorflow.org/guide/ragged\_tensor#tfraggedtensorfrom\_row\_splits .

Let's say we have 5 rays with 2, 1, 3, 0, 2 intersections then the output could look like this

ray_splits : [0, 2, 3, 6, 6, 8] # note that the length of this is num_rays+1
t_hit: [t1, t2, t3, t4, t5, t6, t7, t8]

# ray_splits can be used to iterate over all intersections for each ray
for ray_id, (start, end) in enumerate(zip(ray_splits[:-1], ray_splits[1:])):
    for i,t in enumerate(t_hit[start:end]):
        print(f'ray {ray_id}, intersection {i} at {t}')

We can also think about having both ray_splits and ray_ids to allow iterating over rays and intersections easily.

Copy link
Contributor Author

@dmitrishastin dmitrishastin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks! Didn't realise to look in python folder.

I have added the unit tests; however, while running them, I discovered that the function isn't nearly as fast as I'd hoped - it gets pretty slow in my setting when >250K rays are used. In fact, this increase in running time appears to be exponential:

import numpy as np
import open3d as o3d
import time

cube = o3d.t.geometry.TriangleMesh.from_legacy(
  o3d.geometry.TriangleMesh.create_box())

scene = o3d.t.geometry.RaycastingScene()
scene.add_triangles(cube)

rs = np.random.RandomState(123)

def time_fn(nrays):
  rays = o3d.core.Tensor.from_numpy(rs.rand(nrays, 6).astype(np.float32))
  start = time.time()
  _ = scene.list_intersections(rays)
  finish = time.time()
  return finish - start

rpt_times = 5
nrays = np.array([100, 500, 1000, 5000, 10000, 50000, 100000, 150000, 200000, 250000])
runtime = np.empty((nrays.size, rpt_times))

for r in range(rpt_times):
  for i, n in enumerate(nrays):
    runtime[i, r] = time_fn(n)
    
print(np.round(runtime, 1))
[[ 0.   0.   0.   0.   0. ]
 [ 0.   0.   0.   0.   0. ]
 [ 0.   0.   0.   0.   0. ]
 [ 0.   0.   0.   0.   0. ]
 [ 0.   0.   0.   0.   0. ]
 [ 0.1  0.2  0.2  0.2  0.2]
 [ 0.7  1.   1.   1.   1. ]
 [ 2.9  2.8  2.8  2.7  2.5]
 [ 8.   7.1  8.2  8.3  5.7]
 [16.1 14.1 13.3 12.2 10.8]]

rel_runtime = runtime / nrays[:, None]

print(np.round(rel_runtime, 5))
[[0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [0.e+00 0.e+00 0.e+00 0.e+00 0.e+00]
 [1.e-05 1.e-05 1.e-05 1.e-05 1.e-05]
 [2.e-05 2.e-05 2.e-05 2.e-05 2.e-05]
 [4.e-05 4.e-05 4.e-05 4.e-05 3.e-05]
 [6.e-05 6.e-05 5.e-05 5.e-05 4.e-05]]
 

For now, I have simply set the number of rays in unit tests to be one order less - but will be grateful on any thoughts.

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @benjaminum)


cpp/open3d/t/geometry/RaycastingScene.cpp line 553 at r1 (raw file):

Previously, benjaminum (Benjamin Ummenhofer) wrote…

The & is missing to avoid the copy -> const Eigen::VectorXi& cumsum

Thanks - still getting to grips with c++!


cpp/pybind/t/geometry/raycasting_scene.cpp line 203 at r1 (raw file):

Previously, benjaminum (Benjamin Ummenhofer) wrote…

The idea is to use the cumulative sum for pointing to the start and end of the intersections for each ray similar to this example https://www.tensorflow.org/guide/ragged\_tensor#tfraggedtensorfrom\_row\_splits .

Let's say we have 5 rays with 2, 1, 3, 0, 2 intersections then the output could look like this

ray_splits : [0, 2, 3, 6, 6, 8] # note that the length of this is num_rays+1
t_hit: [t1, t2, t3, t4, t5, t6, t7, t8]

# ray_splits can be used to iterate over all intersections for each ray
for ray_id, (start, end) in enumerate(zip(ray_splits[:-1], ray_splits[1:])):
    for i,t in enumerate(t_hit[start:end]):
        print(f'ray {ray_id}, intersection {i} at {t}')

We can also think about having both ray_splits and ray_ids to allow iterating over rays and intersections easily.

Thanks for the clarification! I have added ray_splits as suggested, and included your example in python documentation. I have left ray_ids in for now - but can take it out if you think that's best.

@benjaminum
Copy link
Contributor

@dmitrishastin I made some minor changes to the example. I like the intersection computation with the primitive_uvs. It is useful code also for interpolating vertex colors and other attributes and I don't think we have such an example yet.

About the shape of the rays, I tested {10,10,6} but the results don't look correct. There seems to be uninitialized memory in the output.
I think it makes sense to drop support for arbitrary dimensions for this function and just allow {n,6} since it is not possible to organize variable length output with multiple dims. Can you have another look at this? Checking for the number of dims and updating the unit test for the shape.

Copy link
Contributor Author

@dmitrishastin dmitrishastin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great - thank you. Using t_hit to compute intersection coordinates is a better method and should have definitely been demonstrated! As you suggest, I am using the code for interpolating surface normals at intersection points for my application, hence the original example.

Regarding shapes, I agree - I was not careful in how I handled arbitrary dimentions.
I have corrected the code - it still supports arbitrary dimentions for input but returns linear outputs as intended.
I have also made this clearer in documentation. The unit test for the shape has been updated.

Quick sanity check:

import open3d as o3d
import numpy as np

scene = o3d.t.geometry.RaycastingScene()
d = o3d.data.MonkeyModel()
_ = scene.add_triangles(o3d.t.io.read_triangle_mesh(d.path))
    
rs = np.random.RandomState(123)
rays = o3d.core.Tensor.from_numpy(rs.rand(4, 5, 6).astype(np.float32))
lx = scene.list_intersections(rays)

print(rays.shape)
print(lx['ray_splits'].shape)
print(lx['t_hit'].shape)
print(lx['ray_splits'])
print(lx['t_hit'])

Output:

SizeVector[4, 5, 6]
SizeVector[21]
SizeVector[12]
[0 0 1 1 2 3 4 5 5 5 6 7 7 8 8 9 9 9 10 11 11]
[0.122634 0.438293 0.0261266 0.305563 0.273778 0.0826897 0.264928 0.0486981 0.0776304 0.348788 0.550925 0.191991]

Reshape:

rays = rays.reshape((20, 6))
lx = scene.list_intersections(rays)
print(rays.shape)
print(lx['ray_splits'].shape)
print(lx['t_hit'].shape)
print(lx['ray_splits'])
print(lx['t_hit'])

Output:

SizeVector[20, 6]
SizeVector[21]
SizeVector[12]
[0 0 1 1 2 3 4 5 5 5 6 7 7 8 8 9 9 9 10 11 11]
[0.122634 0.438293 0.0261266 0.305563 0.273778 0.0826897 0.264928 0.0486981 0.0776304 0.348788 0.550925 0.191991]

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @benjaminum)

Copy link
Contributor

@benjaminum benjaminum left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think only one more change to the documentation is needed and then this PR is ready for merging :)

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @dmitrishastin)


cpp/open3d/t/geometry/RaycastingScene.h line 107 at r6 (raw file):

                                    const int nthreads = 0);

    /// \brief Computes the closest points on the surfaces of the scene.

The documentation ListIntersections and ComputeClosestPoints are mixed up


cpp/open3d/t/geometry/RaycastingScene.h line 128 at r6 (raw file):

            const core::Tensor &rays, const int nthreads = 0);

    /// \brief Lists the intersections of the rays with the scene

The documentation ListIntersections and ComputeClosestPoints are mixed up

Copy link
Contributor Author

@dmitrishastin dmitrishastin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So this is confusing. I believe CountIntersections docs were already incorrect in the main branch.

While preparing my previous commit, I noticed only a part of this and I updated accordingly. After pushing the commit, I realised that the whole CountIntersections text was incorrect and not just that part. I then decided to simply revert documentation to its original text and not touch CountIntersections at all as this had nothing to do with my own PR. Admittedly, this required another commit that I pushed after replying to your last comment, and you may not have seen it.

Overall, I only changed RaycastingScene.h by adding ListIntersections, as is seen in the files section of the PR. I believe that the final documentation for ListIntersection is correct. Please let me know if I've missed anything. I am also happy to correct documentation for CountIntersections as part of my PR - should you deem this appropriate, of course!

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @benjaminum)

Copy link
Contributor

@benjaminum benjaminum left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Documentation is fixed now.
:lgtm:

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @dmitrishastin)

Copy link
Contributor Author

@dmitrishastin dmitrishastin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Many thanks for all your help!
:lgtm:

Reviewable status: 0 of 4 files reviewed, 2 unresolved discussions (waiting on @benjaminum)

@ssheorey ssheorey merged commit 7c0acac into isl-org:master Nov 3, 2023
36 of 37 checks passed
@dmitrishastin dmitrishastin deleted the list_intersections branch November 4, 2023 10:56
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
3 participants