Forum

Author Topic: Why is my custum projection giving different results than camera.project()?  (Read 10837 times)

Abbsalehi

  • Newbie
  • *
  • Posts: 5
    • View Profile
Hi all,

I'm trying to manually implement a projection function in Metashape (python API version 2.2 ) to convert a 3D point (in chunk/world coordinates) into image pixel coordinates, mimicking what camera.project() does internally.

However, I'm getting slightly different results between the built-in camera.project(point) function and my custom implementation. Here's the code I wrote:

Code: [Select]
def manual_project(point_world, camera):
    T = camera.transform.inv()
    point_cam = T.mulp(point_world)

    X, Y, Z = point_cam.x, point_cam.y, point_cam.z

    # Step 1: Normalize
    x = X / Z
    y = Y / Z
    r2 = x**2 + y**2
    r4 = r2**2
    r6 = r2 * r4
    r8 = r4**2

    # Step 2: Get calibration parameters
    calib = camera.sensor.calibration
    f = calib.f
    cx = calib.cx
    cy = calib.cy
    K1, K2, K3, K4 = calib.k1, calib.k2, calib.k3, calib.k4
    P1, P2 = calib.p1, calib.p2
    B1, B2 = calib.b1, calib.b2
    pixel_size_x = camera.sensor.pixel_size[0]
    pixel_size_y = camera.sensor.pixel_size[1]
    w, h = camera.sensor.width, camera.sensor.height

    # Step 3: Apply distortion
    radial = 1 + K1*r2 + K2*r4 + K3*r6 + K4*r8
    x_dist = x * radial + (P1 * (r2 + 2*x**2) + 2 * P2 * x * y)
    y_dist = y * radial + (P2 * (r2 + 2*y**2) + 2 * P1 * x * y)

    # Step 4: Convert to pixel coordinates
    u = w * 0.5 + cx + x_dist * f + x_dist * B1 + y_dist * B2
    v = h * 0.5 + cy + y_dist * f

    return np.array([u, v])

point_world = Vector([-22.261804864818998, 5.591732564987695, 10.145055367315823]) #3D Point in Chunk Coordinates

f  = 3945.7686579791701
cx = 10.60110469705829
cy = -60.54494592517153

K1 = -0.09472477721894505
K2 =  0.09189283628424988
K3 = -0.023947156862946
K4 =  0.0

P1 =  0.0004860575870856356
P2 = -0.00018801209260916065

B1 = 0.0
B2 = 0.0

image size = 5760 x 3840


Is there something wrong with how I'm applying here? please let me know if you need anything else thanks
Any insights on aligning the behavior of my function with camera.project() would be greatly appreciated!

Thanks in advance 🙏
« Last Edit: March 28, 2025, 03:48:32 PM by Abbsalehi »

Paulo

  • Hero Member
  • *****
  • Posts: 1599
    • View Profile
Hello Abbsalehi,

The world to camera coordinates part is not correct..

Given crs = chunk.crs, T = chunk.transform.matrix and CT = camera.transform then

point_cam = CT.inv().mulp(T.inv().mulp(crs.unproject(point_world))).


PS. if by world coordinates you mean the internal chunk coordinates then your code should be correct...

i.e. point_cam = CT.inv().mulp(point_int)
« Last Edit: March 27, 2025, 08:29:53 PM by Paulo »
Best Regards,
Paul Pelletier,
Surveyor

Abbsalehi

  • Newbie
  • *
  • Posts: 5
    • View Profile
Hello Abbsalehi,

The world to camera coordinates part is not correct..

Given crs = chunk.crs, T = chunk.transform.matrix and CT = camera.transform then

point_cam = CT.inv().mulp(T.inv().mulp(crs.unproject(point_world))).


PS. if by world coordinates you mean the internal chunk coordinates then your code should be correct...

i.e. point_cam = CT.inv().mulp(point_int)

Hi Paulo, and thanks for your reply; I am importing my point cloud like this below code; is it in the internal chunk coordinates?

Code: [Select]
doc_2.open("metshape_file.psx")
chunk = doc_2.chunk
chunk = doc_2.chunk
point_cloud = chunk.point_cloud
« Last Edit: March 27, 2025, 09:56:02 PM by Abbsalehi »

Paulo

  • Hero Member
  • *****
  • Posts: 1599
    • View Profile
Hi Abbs,

yes the coordinates of point cloud are in internal ....
Best Regards,
Paul Pelletier,
Surveyor

Abbsalehi

  • Newbie
  • *
  • Posts: 5
    • View Profile
Hi Abbs,

yes the coordinates of point cloud are in internal ....

I tested your line of code. Also, and it is still off, I think something is off here; I used the frame type camera equations from the appendix section https://www.agisoft.com/pdf/metashape-pro_2_2_en.pdf. Do you have any other thoughts?
« Last Edit: March 27, 2025, 10:46:39 PM by Abbsalehi »

Paulo

  • Hero Member
  • *****
  • Posts: 1599
    • View Profile
If you use internal coordinates then your original code works... see following:

Code: [Select]
c,pt
Out[16]: 2025-03-27 15:05:03
2025-03-27 15:05:03 (<Camera 'frame_315.jpg'>,
2025-03-27 15:05:03  Vector([-11.877580631661171, -43.58721543611359, 1.7677075429534725]))

manual_project(pt,c)
Out[17]: 2025-03-27 15:05:09 array([327.84172962, 228.84602583])

c.project(pt)
Out[18]: 2025-03-27 15:05:17 Vector([327.8417296236855, 228.84602582719452])
Best Regards,
Paul Pelletier,
Surveyor

Abbsalehi

  • Newbie
  • *
  • Posts: 5
    • View Profile
If you use internal coordinates then your original code works... see following:

Code: [Select]
c,pt
Out[16]: 2025-03-27 15:05:03
2025-03-27 15:05:03 (<Camera 'frame_315.jpg'>,
2025-03-27 15:05:03  Vector([-11.877580631661171, -43.58721543611359, 1.7677075429534725]))

manual_project(pt,c)
Out[17]: 2025-03-27 15:05:09 array([327.84172962, 228.84602583])

c.project(pt)
Out[18]: 2025-03-27 15:05:17 Vector([327.8417296236855, 228.84602582719452])

Thanks for your help!

Abbsalehi

  • Newbie
  • *
  • Posts: 5
    • View Profile
If you use internal coordinates then your original code works... see following:

Code: [Select]
c,pt
Out[16]: 2025-03-27 15:05:03
2025-03-27 15:05:03 (<Camera 'frame_315.jpg'>,
2025-03-27 15:05:03  Vector([-11.877580631661171, -43.58721543611359, 1.7677075429534725]))

manual_project(pt,c)
Out[17]: 2025-03-27 15:05:09 array([327.84172962, 228.84602583])

c.project(pt)
Out[18]: 2025-03-27 15:05:17 Vector([327.8417296236855, 228.84602582719452])

Thanks for your help!

Hi Paulo,

I just projected the `points_int` and overlaid the projected points on the image, and I received the result shown in the screenshot below. The highlighted area does not match what I obtained using `camera.project`. I'm not sure what might be wrong with my code. Do you have any thoughts on this?

Thanks!

Paulo

  • Hero Member
  • *****
  • Posts: 1599
    • View Profile
Both methods will give same result however you should check visibility between 3d internal point and respective camera...

For more details, please send me a PM
Best Regards,
Paul Pelletier,
Surveyor