Projecting from 3D point coordinates to the best areal camera/image, save that.

[Emanuele1234]

Hello everyone,

I start a project where I need to project some 3D coordinates ( that I don't have on metashape, but in a csv.file) on 2d areal camera/image using python linked to a metashape project. Then I need to cut off and save the best areal image were each point is projected, so I will have a folder with all the best image of the single point in the map.

Someone can help me?

I started using this script that I saw in another post here in the forum:

Code: [Select]

import Metashape, math, csv

path = "I:/Backsjon190527/save_progress/"
doc = Metashape.Document()
doc.open(path + "Finale.psx")
chunk = doc.chunk
crs = chunk.crs
print(crs)

# export format:
# point_ID X Y Z
# [camera_label x_proj y_proj]

M = chunk.transform.matrix
point_cloud = chunk.point_cloud
projections = point_cloud.projections
points = point_cloud.points
npoints = len(points)
print(npoints)
tracks = point_cloud.tracks



path1 = Metashape.app.getSaveFileName(path + 'New.csv',
                                     filter="Text / CSV (*.txt *.csv);;All files (*.*)")
file = open("New.csv", 'wt')

print("Script started...")

point_ids = [-1] * len(point_cloud.tracks)
for point_id in range(0, npoints):
    point_ids[points[point_id].track_id] = point_id
points_proj = {}

for photo in chunk.cameras:

    if not photo.transform:
        continue
    T = photo.transform.inv()
    calib = photo.sensor.calibration

    for proj in projections[photo]:
        track_id = proj.track_id
        point_id = point_ids[track_id]

        if point_id < 0:
            continue
        if not points[point_id].valid:  # skipping invalid points
            continue

        point_index = point_id
        if point_index in points_proj.keys():
            x, y = proj.coord
            points_proj[point_index] = (points_proj[point_index] + "\n" + photo.label + "\t{:.2f}\t{:.2f}".format(x, y))

        else:
            x, y = proj.coord
            points_proj[point_index] = ("\n" + photo.label + "\t{:.2f}\t{:.2f}".format(x, y))

for point_index in range(npoints):

    if not points[point_index].valid:
        continue

    coord = M * points[point_index].coord
    coord.size = 3
    if chunk.crs:
        # coord
        X, Y, Z = chunk.crs.project(coord)
    else:
        X, Y, Z = coord

    line = points_proj[point_index]

    file.write("{}\t{:.6f}\t{:.6f}\t{:.6f}\t{:s}\n".format(point_index, X, Y, Z, line))


file.close()
print("Finished")
This code probably is not good for what I want to do...

Thank you in advance for who will help me!
Best regard.
Emanuele.

[Paulo]

Hi Emanuele,

suppose you have a point P whose coordinates (X,Y,Z) are in same CRS as chunk, then following code will populate a dictionary cam_dist_P with distance from projected point to camera center or Principal point for each camera where it appears:

Code: [Select]

chunk = Metashape.app.document.chunk
T = chunk.transform.matrix

selcam = list()
mindis = 1e09
P = Metashape.Vector((459307.345426, 7192625.624622, 158.808)) # vector with X,Y,Z coordinates of point in same CRS as chunk
cam_dist_P = {}    # dict containing  distance in pixels of P to each camera Principal Point

for camera in chunk.cameras:
sensor = camera.sensor
w = sensor.width
h = sensor.height
cx = sensor.calibration.cx  # Principal point coordinates x
cy = sensor.calibration.cy  # Principal point coordinates y
if not camera.project(T.inv().mulp(chunk.crs.unproject(P))): # check if projected Point within sensor size
continue
xp = camera.project(T.inv().mulp(chunk.crs.unproject(P))).x  # Prrojected point coordinates x
yp = camera.project(T.inv().mulp(chunk.crs.unproject(P))).y  # Prrojected point coordinates y
if xp < 0 or xp > w or yp < 0 or yp > h: # check if projected PointI within sensor size
continue
dist = (camera.project(T.inv().mulp(chunk.crs.unproject(P)))-Metashape.Vector((w/2+cx,h/2+cy))).norm() # distance in pixels from projected point and PP
cam_dist_P[camera.label] = dist
if dist < mindis :
mindis = dist
mincam = camera # camera with min distance to Principal Point
mincam.selected = True
selcam.append(mincam)
so cam_dist_ P would contain:

Code: [Select]

{'DJI_0357.tiff': 82.53159973111778,
 'DJI_0359.tiff': 88.38535725791587,
 'DJI_0361.tiff': 123.93762482172889,
 'DJI_0363.tiff': 130.49876994869547,
 'DJI_0365.tiff': 121.34328558462049,
 'DJI_0367.tiff': 188.1538982469012,
 'DJI_0369.tiff': 245.78868399774544}
and sel_cam would contain camera with minimum distance to projected point  [<Camera 'DJI_0357.tiff'>]

Hope this can orient your pursuit....

[Emanuele1234]

Hi Paulo,
thank you for answer me!

Probably the script can help me! The only problem is that before I find the best image per marker, I need to project the 3D coordinates for that point to 2D coordinates. So I can use 2D coordinates to find the best image and cut off a buffer of that around the marker. But the script you give me, need vector with 3 coordinates to work. Do you think there is a method to use 2D coordinates point? Or in alternative, Do you think can be better try to use 3D coordinates to find the best Image and then transform that in 2D coordinates?... What do you think?

Thank you again for your help!
Best regards,
Emanuele.

[Paulo]

Emanuele,

this is exactly what it does. Given a point P with 3d coordinates (X,Y,Z) in chunk CRS, then the projected image pixel coordinates for a given camera is given by:

xp, yp = camera.project(T.inv().mulp(chunk.crs.unproject(Metashape.Vector((X,Y,Z))))).

So in the script, if you want the best camera projected coordinates (one with projection closest to camera center) then they would be given by:

xp, yp = mincam.project(T.inv().mulp(chunk.crs.unproject(Metashape.Vector((X,Y,Z))))).

[Emanuele1234]

Ah, oh my god, thank you so mutch! I'm sorry do not undestand before! I'm just a beginner in python, so a lot of things are really difficult for me! But thank you very mutch!

[Paulo]

No worry Emanuele, It took quite a while for me to get a grip on different coordinates sytems used in MS and their transformations!