rnd-team-dev / plotoptix

Data visualisation and ray tracing in Python based on OptiX 7.7 framework.
https://rnd.team/plotoptix
Other
499 stars 26 forks source link

How can I customize the ray into the program and calculate intersection #43

Closed mshmoon closed 1 year ago

mshmoon commented 1 year ago

How can I customize the ray into the program and calculate intersection

mshmoon commented 1 year ago

To be specific, I have many source points, and each source point has many directions. I want to transfer the rays to the GPU for computing intersection.

robertsulej commented 1 year ago

You can follow this discussion and linked examples. There are also related questions in issues 36, 38, 39, and 41 that we were discussing recently.

Have a look at code samples in these issues, ask if you need help :)

mshmoon commented 1 year ago

You can follow this discussion and linked examples. There are also related questions in issues 36, 38, 39, and 41 that we were discussing recently.

Have a look at code samples in these issues, ask if you need help :)

Thank you very much for your reply. I want to get some answers as soon as possible, so I will express my questions more specifically.

If: I have a these ray origin [[0,0,0], [1,1,1], [2,2,2] The direction of the ray is [[3,3,3], [4,4,4], [5,5,5]].

I want to know if I can put the origin and direction into the program to calculate the intersection. My current idea is to set two 1-dimensional textures, then [[0,0,0], [1,1,1], [2,2,2]] as texture A, and [[3,3,3], [4,4,4], [5,5,5]] as texture B, and send them into the program for intersection

robertsulej commented 1 year ago

This is exactly what this sample is doing.

Things to note (all is in the sample, just to remember):

mshmoon commented 1 year ago

This is exactly what this sample is doing.

Things to note (all is in the sample, just to remember):

  • origins and directons are in 2D textures,
  • rendering size need to be same as textures above
  • do not forget to wait for gpu.

Thank you, if I have 10000000 rays, it will be work?

mshmoon commented 1 year ago

This is exactly what this sample is doing.

Things to note (all is in the sample, just to remember):

  • origins and directons are in 2D textures,
  • rendering size need to be same as textures above
  • do not forget to wait for gpu.

After your guidance, my program seems to be running. Thank you very much. But I have a question. The distance calculated by this program seems to have some errors. For example, the distance should be an integer 1, but the resulting distance is 0.998. How can we avoid this problem

mshmoon commented 1 year ago

@robertsulej
It seems that there are some problems with ray sampling, such as:

this is my origin textures: [[[ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.]]

[[ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.] [ 5. 25. 125. 0.]]]

this is my directions textures [[[ 0.66727098 -0.22638111 -0.70957807 0. ] [ 0.64108496 -0.23752188 -0.72978999 0. ] [ 0.61411788 -0.24837328 -0.74911277 0. ] [ 0.58640259 -0.25892207 -0.76752287 0. ] [ 0.55797287 -0.2691554 -0.78499787 0. ] [ 0.52886333 -0.27906081 -0.80151646 0. ]]

[[ 0.49910946 -0.28862622 -0.81705853 0. ] [ 0.46874751 -0.29783999 -0.83160514 0. ] [ 0.43781445 -0.30669089 -0.84513857 0. ] [ 0.40634799 -0.31516813 -0.85764233 0. ] [ 0.37438645 -0.32326139 -0.86910118 0. ] [ 0.34196878 -0.3309608 -0.87950117 0. ]]],

From the above data, it is impossible to sample the same rays from these two matrices,but two identical rays are sampled,Please see the following two key lines (the data is come from rt._hit_pos: x,y,z,distant format ): x y z distant 39.998 9.5610695 79.22586 59.6531 39.998 4.7809324 67.75382 70.077065 39.998 12.032505 85.15715 54.593666 39.998 8.116671 75.759384 62.725952 39.998 9.5610695 79.22586 59.6531 39.998 4.7809324 67.75382 70.077065 39.998 6.5132847 71.911354 66.219505 4353428800000000.0 -3052965300000000.0 -8411762700000000.0 9951416000000000.0 4687475000000000.0 -2978400000000000.0 -8316051500000000.0 1e+16 3743864700000000.0 -3232613800000000.0 -8691012000000000.0 1e+16 39.998 6.5132847 71.911354 66.219505 4353428800000000.0 -3052965300000000.0 -8411762700000000.0 9951416000000000.0

Thanks for your help

mshmoon commented 1 year ago

This is exactly what this sample is doing.

Things to note (all is in the sample, just to remember):

  • origins and directons are in 2D textures,
  • rendering size need to be same as textures above
  • do not forget to wait for gpu.

I found that the result is more accurate when the texture is set to square texture(for example: 4x4,8x8,16x16). This conclusion surprised me,why?

robertsulej commented 1 year ago

Hi!

As for the number of rays, 10M is not that much, it is just a 3K image. You can go much higher.

I double checked results of distance calculations and they look correct. Eg the code below is shooting parallel rays along Z axis towards a plane at a distance 1.0:

import numpy as np

from plotoptix import NpOptiX
from plotoptix.materials import m_flat

import threading

done = threading.Event()

def wait_for_gpu(rt: NpOptiX) -> None:
    print("gpu done")
    done.set()

def main():
    nu = 17
    nv = 17

    rt = NpOptiX(on_launch_finished=wait_for_gpu, width=nu, height=nv)
    rt.set_param(min_accumulation_step=1, max_accumulation_frames=1)
    rt.setup_material("flat", m_flat)

    rt.set_data("surf", pos=[-1.5,-1.5,0], u=[0,3,0], v=[3,0,0], geom="Parallelograms", mat="flat")

    u = np.linspace(-2, 2, nu)
    v = np.linspace(2, -2, nv)

    V, U = np.meshgrid(v, u)
    W = np.full((nu, nv), -1)

    src_tex = np.stack((U, V, W, np.zeros((nu, nv)))).T
    rt.set_texture_2d("src", src_tex)

    cx = np.zeros((nu, nv))
    cy = np.zeros((nu, nv))
    cz = np.ones((nu, nv))
    r = np.full((nu, nv), 10)
    dir_tex = np.stack((cx, cy, cz, r)).T
    rt.set_texture_2d("dir", dir_tex)

    rt.setup_camera("custom_cam", cam_type="CustomProjXYZtoDir", textures=["src", "dir"])

    rt.start()

    print("working...")

    if done.wait(10):
        f = rt._geo_id[:,:,1].reshape(rt._height, rt._width)
        d = rt._hit_pos[:,:,3].reshape(rt._height, rt._width)
        d[f == 0xFFFFFFFF] = 0

        print(d)

    else:
        print("timeout")

    rt.close()

if __name__ == '__main__':
    main()

You should see a trace of the flat face:

[[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. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 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. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]]

A bit more complex code is shooting a cone of rays from a single point towards a plane:

import numpy as np

from plotoptix import NpOptiX
from plotoptix.materials import m_flat

import threading

done = threading.Event()

def wait_for_gpu(rt: NpOptiX) -> None:
    print("gpu done")
    done.set()

def main():
    np.set_printoptions(precision=3, linewidth=150)

    nu = 17
    nv = 17

    rt = NpOptiX(on_launch_finished=wait_for_gpu, width=nu, height=nv)
    rt.set_param(min_accumulation_step=1, max_accumulation_frames=1)
    rt.setup_material("flat", m_flat)

    rt.set_data("surf", pos=[-3,-3,0], u=[0,6,0], v=[6,0,0], geom="Parallelograms", mat="flat")

    src = np.zeros((nu, nv, 3))
    src[:,:,2] = -1
    src_tex = np.c_[src, np.zeros((nu, nv, 1))]
    rt.set_texture_2d("src", src_tex)

    u = np.linspace(-2, 2, nu)
    v = np.linspace(2, -2, nv)

    V, U = np.meshgrid(v, u)
    W = np.zeros((nu, nv))

    dst = np.stack((U, V, W)).T

    dir_tex = dst - src
    dir_tex = dir_tex / np.linalg.norm(dir_tex, axis=-1)[:, :, np.newaxis]
    print(dir_tex)
    dir_tex = np.c_[dir_tex, np.full((nu, nv, 1), 100)]
    rt.set_texture_2d("dir", dir_tex)

    rt.setup_camera("custom_cam", cam_type="CustomProjXYZtoDir", textures=["src", "dir"])

    rt.start()

    print("working...")

    if done.wait(10):
        f = rt._geo_id[:,:,1].reshape(rt._height, rt._width)
        d = rt._hit_pos[:,:,3].reshape(rt._height, rt._width)
        d[f == 0xFFFFFFFF] = 0

        print(d)
        print(d[8,8])

    else:
        print("timeout")

    rt.close()

if __name__ == '__main__':
    main()

And also here all the distances are exactly as expected (1.0 in the center, 3.0 in the corners):

[[3.    2.839 2.693 2.562 2.449 2.358 2.291 2.25  2.236 2.25  2.291 2.358 2.449 2.562 2.693 2.839 3.   ]
 [2.839 2.669 2.512 2.372 2.25  2.151 2.077 2.031 2.016 2.031 2.077 2.151 2.25  2.372 2.512 2.669 2.839]
 [2.693 2.512 2.345 2.194 2.062 1.953 1.871 1.82  1.803 1.82  1.871 1.953 2.062 2.194 2.345 2.512 2.693]
 [2.562 2.372 2.194 2.031 1.887 1.768 1.677 1.62  1.601 1.62  1.677 1.768 1.887 2.031 2.194 2.372 2.562]
 [2.449 2.25  2.062 1.887 1.732 1.601 1.5   1.436 1.414 1.436 1.5   1.601 1.732 1.887 2.062 2.25  2.449]
 [2.358 2.151 1.953 1.768 1.601 1.458 1.346 1.275 1.25  1.275 1.346 1.458 1.601 1.768 1.953 2.151 2.358]
 [2.291 2.077 1.871 1.677 1.5   1.346 1.225 1.146 1.118 1.146 1.225 1.346 1.5   1.677 1.871 2.077 2.291]
 [2.25  2.031 1.82  1.62  1.436 1.275 1.146 1.061 1.031 1.061 1.146 1.275 1.436 1.62  1.82  2.031 2.25 ]
 [2.236 2.016 1.803 1.601 1.414 1.25  1.118 1.031 1.    1.031 1.118 1.25  1.414 1.601 1.803 2.016 2.236]
 [2.25  2.031 1.82  1.62  1.436 1.275 1.146 1.061 1.031 1.061 1.146 1.275 1.436 1.62  1.82  2.031 2.25 ]
 [2.291 2.077 1.871 1.677 1.5   1.346 1.225 1.146 1.118 1.146 1.225 1.346 1.5   1.677 1.871 2.077 2.291]
 [2.358 2.151 1.953 1.768 1.601 1.458 1.346 1.275 1.25  1.275 1.346 1.458 1.601 1.768 1.953 2.151 2.358]
 [2.449 2.25  2.062 1.887 1.732 1.601 1.5   1.436 1.414 1.436 1.5   1.601 1.732 1.887 2.062 2.25  2.449]
 [2.562 2.372 2.194 2.031 1.887 1.768 1.677 1.62  1.601 1.62  1.677 1.768 1.887 2.031 2.194 2.372 2.562]
 [2.693 2.512 2.345 2.194 2.062 1.953 1.871 1.82  1.803 1.82  1.871 1.953 2.062 2.194 2.345 2.512 2.693]
 [2.839 2.669 2.512 2.372 2.25  2.151 2.077 2.031 2.016 2.031 2.077 2.151 2.25  2.372 2.512 2.669 2.839]
 [3.    2.839 2.693 2.562 2.449 2.358 2.291 2.25  2.236 2.25  2.291 2.358 2.449 2.562 2.693 2.839 3.   ]]

If you set the size of the plane matching the cone of rays (rt.set_data("surf", pos=[-2,-2,0], u=[0,4,0], v=[4,0,0], geom="Parallelograms", mat="flat")), you can expect single rays missing the plane on its edges, but this is normal effect of the precision.

If you share your code, I can have a look (sorry for delay, I am replying on the first free time slot I have).

mshmoon commented 1 year ago

Hi!

As for the number of rays, 10M is not that much, it is just a 3K image. You can go much higher.

I double checked results of distance calculations and they look correct. Eg the code below is shooting parallel rays along Z axis towards a plane at a distance 1.0:

import numpy as np

from plotoptix import NpOptiX
from plotoptix.materials import m_flat

import threading

done = threading.Event()

def wait_for_gpu(rt: NpOptiX) -> None:
    print("gpu done")
    done.set()

def main():
    nu = 17
    nv = 17

    rt = NpOptiX(on_launch_finished=wait_for_gpu, width=nu, height=nv)
    rt.set_param(min_accumulation_step=1, max_accumulation_frames=1)
    rt.setup_material("flat", m_flat)

    rt.set_data("surf", pos=[-1.5,-1.5,0], u=[0,3,0], v=[3,0,0], geom="Parallelograms", mat="flat")

    u = np.linspace(-2, 2, nu)
    v = np.linspace(2, -2, nv)

    V, U = np.meshgrid(v, u)
    W = np.full((nu, nv), -1)

    src_tex = np.stack((U, V, W, np.zeros((nu, nv)))).T
    rt.set_texture_2d("src", src_tex)

    cx = np.zeros((nu, nv))
    cy = np.zeros((nu, nv))
    cz = np.ones((nu, nv))
    r = np.full((nu, nv), 10)
    dir_tex = np.stack((cx, cy, cz, r)).T
    rt.set_texture_2d("dir", dir_tex)

    rt.setup_camera("custom_cam", cam_type="CustomProjXYZtoDir", textures=["src", "dir"])

    rt.start()

    print("working...")

    if done.wait(10):
        f = rt._geo_id[:,:,1].reshape(rt._height, rt._width)
        d = rt._hit_pos[:,:,3].reshape(rt._height, rt._width)
        d[f == 0xFFFFFFFF] = 0

        print(d)

    else:
        print("timeout")

    rt.close()

if __name__ == '__main__':
    main()

You should see a trace of the flat face:

[[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. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 0. 0.]
 [0. 0. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 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. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]]

A bit more complex code is shooting a cone of rays from a single point towards a plane:

import numpy as np

from plotoptix import NpOptiX
from plotoptix.materials import m_flat

import threading

done = threading.Event()

def wait_for_gpu(rt: NpOptiX) -> None:
    print("gpu done")
    done.set()

def main():
    np.set_printoptions(precision=3, linewidth=150)

    nu = 17
    nv = 17

    rt = NpOptiX(on_launch_finished=wait_for_gpu, width=nu, height=nv)
    rt.set_param(min_accumulation_step=1, max_accumulation_frames=1)
    rt.setup_material("flat", m_flat)

    rt.set_data("surf", pos=[-3,-3,0], u=[0,6,0], v=[6,0,0], geom="Parallelograms", mat="flat")

    src = np.zeros((nu, nv, 3))
    src[:,:,2] = -1
    src_tex = np.c_[src, np.zeros((nu, nv, 1))]
    rt.set_texture_2d("src", src_tex)

    u = np.linspace(-2, 2, nu)
    v = np.linspace(2, -2, nv)

    V, U = np.meshgrid(v, u)
    W = np.zeros((nu, nv))

    dst = np.stack((U, V, W)).T

    dir_tex = dst - src
    dir_tex = dir_tex / np.linalg.norm(dir_tex, axis=-1)[:, :, np.newaxis]
    print(dir_tex)
    dir_tex = np.c_[dir_tex, np.full((nu, nv, 1), 100)]
    rt.set_texture_2d("dir", dir_tex)

    rt.setup_camera("custom_cam", cam_type="CustomProjXYZtoDir", textures=["src", "dir"])

    rt.start()

    print("working...")

    if done.wait(10):
        f = rt._geo_id[:,:,1].reshape(rt._height, rt._width)
        d = rt._hit_pos[:,:,3].reshape(rt._height, rt._width)
        d[f == 0xFFFFFFFF] = 0

        print(d)
        print(d[8,8])

    else:
        print("timeout")

    rt.close()

if __name__ == '__main__':
    main()

And also here all the distances are exactly as expected (1.0 in the center, 3.0 in the corners):

[[3.    2.839 2.693 2.562 2.449 2.358 2.291 2.25  2.236 2.25  2.291 2.358 2.449 2.562 2.693 2.839 3.   ]
 [2.839 2.669 2.512 2.372 2.25  2.151 2.077 2.031 2.016 2.031 2.077 2.151 2.25  2.372 2.512 2.669 2.839]
 [2.693 2.512 2.345 2.194 2.062 1.953 1.871 1.82  1.803 1.82  1.871 1.953 2.062 2.194 2.345 2.512 2.693]
 [2.562 2.372 2.194 2.031 1.887 1.768 1.677 1.62  1.601 1.62  1.677 1.768 1.887 2.031 2.194 2.372 2.562]
 [2.449 2.25  2.062 1.887 1.732 1.601 1.5   1.436 1.414 1.436 1.5   1.601 1.732 1.887 2.062 2.25  2.449]
 [2.358 2.151 1.953 1.768 1.601 1.458 1.346 1.275 1.25  1.275 1.346 1.458 1.601 1.768 1.953 2.151 2.358]
 [2.291 2.077 1.871 1.677 1.5   1.346 1.225 1.146 1.118 1.146 1.225 1.346 1.5   1.677 1.871 2.077 2.291]
 [2.25  2.031 1.82  1.62  1.436 1.275 1.146 1.061 1.031 1.061 1.146 1.275 1.436 1.62  1.82  2.031 2.25 ]
 [2.236 2.016 1.803 1.601 1.414 1.25  1.118 1.031 1.    1.031 1.118 1.25  1.414 1.601 1.803 2.016 2.236]
 [2.25  2.031 1.82  1.62  1.436 1.275 1.146 1.061 1.031 1.061 1.146 1.275 1.436 1.62  1.82  2.031 2.25 ]
 [2.291 2.077 1.871 1.677 1.5   1.346 1.225 1.146 1.118 1.146 1.225 1.346 1.5   1.677 1.871 2.077 2.291]
 [2.358 2.151 1.953 1.768 1.601 1.458 1.346 1.275 1.25  1.275 1.346 1.458 1.601 1.768 1.953 2.151 2.358]
 [2.449 2.25  2.062 1.887 1.732 1.601 1.5   1.436 1.414 1.436 1.5   1.601 1.732 1.887 2.062 2.25  2.449]
 [2.562 2.372 2.194 2.031 1.887 1.768 1.677 1.62  1.601 1.62  1.677 1.768 1.887 2.031 2.194 2.372 2.562]
 [2.693 2.512 2.345 2.194 2.062 1.953 1.871 1.82  1.803 1.82  1.871 1.953 2.062 2.194 2.345 2.512 2.693]
 [2.839 2.669 2.512 2.372 2.25  2.151 2.077 2.031 2.016 2.031 2.077 2.151 2.25  2.372 2.512 2.669 2.839]
 [3.    2.839 2.693 2.562 2.449 2.358 2.291 2.25  2.236 2.25  2.291 2.358 2.449 2.562 2.693 2.839 3.   ]]

If you set the size of the plane matching the cone of rays (rt.set_data("surf", pos=[-2,-2,0], u=[0,4,0], v=[4,0,0], geom="Parallelograms", mat="flat")), you can expect single rays missing the plane on its edges, but this is normal effect of the precision.

If you share your code, I can have a look (sorry for delay, I am replying on the first free time slot I have).

Thank you for your reply. However, I feel that Plotoptix has little advantage in running speed on cuda, and the computing speed is even slower than open3d on cpu. Maybe there is a problem with my code, so I also want to ask you whether Plotoptix has advantages in computing speed

robertsulej commented 1 year ago

There is a large number of factors affecting speed comparison: number of rays, complexity of the geometry, compilation time of the scene, and of course the python code which can easyli ruin all the performance with a single loop. A lot depends also on what is your aim. Pure calculation of ray intersections with a single plane is trivial and can run very fast on cpu. But this will change when there are thousands-milions of triangles. PlotOptiX is not dedicated primarily to do vary basic calculations alone. Distances, normals, etc are available as a side product of running more complex, usually multi-segment tracing and shading. If you need only distances or hit positions, then it can still be faster than cpu code, but not on 16x16 image and a single plane. There is overhead of the full machinery. It can be set up to some minimal config, but it is there. I can have a look at you code if you like.

mshmoon commented 1 year ago

There is a large number of factors affecting speed comparison: number of rays, complexity of the geometry, compilation time of the scene, and of course the python code which can easyli ruin all the performance with a single loop. A lot depends also on what is your aim. Pure calculation of ray intersections with a single plane is trivial and can run very fast on cpu. But this will change when there are thousands-milions of triangles. PlotOptiX is not dedicated primarily to do vary basic calculations alone. Distances, normals, etc are available as a side product of running more complex, usually multi-segment tracing and shading. If you need only distances or hit positions, then it can still be faster than cpu code, but not on 16x16 image and a single plane. There is overhead of the full machinery. It can be set up to some minimal config, but it is there. I can have a look at you code if you like.

import numpy as np from plotoptix import NpOptiX import matplotlib.pyplot as plt import open3d as o3d import open3d.visualization as vis from building_data import init_env # def least_square(num): for i in range(10000): if i*i>num: return i

def vstack(i,x): need_pad = i * i - np.shape(x)[0] padding = np.zeros((need_pad,3)) x = np.vstack((x,padding)) return x

def hstack(x): need_pad = np.shape(x)[0] padding = np.zeros((need_pad,1)) x = np.hstack((x,padding)) return x

v,f = init_env() # my function build_mesh = o3d.geometry.TriangleMesh() build_mesh.vertices = o3d.utility.Vector3dVector(v) build_mesh.triangles = o3d.utility.Vector3iVector(f)

import time v, t = init_env() print("start load data") rays=np.loadtxt("rays.txt").astype(np.float32) print("load data over")

origins = rays[:,0:3] directions = rays[:,3:6]

num = np.shape(rays)[0] print("rays num:",num) ls_num = least_square(num)

hstack_pad_num = np.zeros((ls_num * ls_num,1))

origins = vstack(ls_num,origins) directions = vstack(ls_num,directions)

origins = hstack(origins).reshape(ls_num ,ls_num ,4) directions = hstack(directions).reshape(ls_num ,ls_num ,4)

s1 = time.time() def displayResults(rt): print("Launch finished.") s2 = time.time() print(s2 - s1) hitPositionsData = rt._hit_pos xHitPositions = hitPositionsData[:, :, 0].reshape(-1) yHitPositions = hitPositionsData[:, :, 1].reshape(-1) zHitPositions = hitPositionsData[:, :, 2].reshape(-1) dHitPositions = hitPositionsData[:, :, 3].reshape(-1)

points= []
for n,d in enumerate(dHitPositions):
    if d<1e+5:
        points.append([xHitPositions[n], yHitPositions[n], zHitPositions[n]])

points = np.array(points)
points_pcd = o3d.geometry.PointCloud()
points_pcd.points = o3d.utility.Vector3dVector(points)
vis.draw([points_pcd,build_mesh])

rt = NpOptiX(on_rt_accum_done=displayResults, width = ls_num , height = ls_num ) rt.set_mesh("Mesh", v,t) rt.set_texture_2d("origins", origins) rt.set_texture_2d("directions", directions)

rt.setup_camera("custom_cam", cam_type="CustomProjXYZtoDir", textures=["origins", "directions"]) rt.set_param(max_accumulation_frames=1) rt.set_uint("path_seg_range", 1, 4) rt.start()

This is my code. Under the same mesh scene and ray number, plotopix is slower than open3d. Is the underlying computing engine of plotoptix nvidia's optix?

robertsulej commented 1 year ago

Sorry for the long delay... if you (or others) are still reading the thread, here are some answers.

Yes, plotoptix is a set of CUDA programs in running OptiX framework.

I could not run exactly your code, since it does not include func generating data. Some quick notes from reading the code: