How to calculate camera y-shift value to compensate when removing vertical perspective distortion

Question

I am trying to write a script to remove the vertical perspective from a camera relative to the center of the view frame (e.g. I want all vertical lines to be parallel when viewed from the camera). In order to achieve this, I set the camera's local x rotation to 90, and then need to compensate using the y-shift is order to keep what was exactly at the center of the frame in the same position.

I am struggling with figuring out how to determine the proper y-shift to achieve this. This is my current attempt, which unfortunately does not work:

def auto_remove_camera_vertical_perspective(cam_obj):
    # ==Add plane and align to cam frustum==
    bpy.ops.mesh.primitive_plane_add(
        enter_editmode=False,
        align="WORLD",
        location=(0, 0, 0),
        scale=(1, 1, 1),
    )
    view_plane = bpy.context.object

    # Align rotation
    bpy.ops.object.constraint_add(type="COPY_ROTATION")
    view_plane.constraints["Copy Rotation"].target = cam_obj
    bpy.ops.constraint.apply(constraint="Copy Rotation", owner="OBJECT")

    # Align location using the cam view frame plane coords
    cam_data = cam_obj.data

    matrix = cam_obj.matrix_world.normalized()
    frame_points = [matrix @ v for v in cam_data.view_frame(scene=bpy.context.scene)]

    x = sum(p.x for p in frame_points) / 4
    y = sum(p.y for p in frame_points) / 4
    z = sum(p.z for p in frame_points) / 4

    center_point = Vector((x, y, z))
    view_plane.location = center_point

    # Scale plane
    scale_factor = 10_000
    view_plane.scale = Vector([scale_factor] * 3)

    # Apply transforms
    bpy.ops.object.transform_apply(location=True, rotation=True, scale=True)

    # ==Ray cast to plane==
    # Get the camera's location and rotation
    ray_origin = cam_obj.location
    cam_rot = cam_obj.rotation_euler.to_matrix()

    # Calculate the ray direction from camera's facing
    ray_direction = cam_rot @ Vector((0, 0, -1))

    # Cast the ray from the camera origin onto the object
    ray_cast_result = view_plane.ray_cast(ray_origin, ray_direction)

    ray_cast_result = view_plane.ray_cast(ray_origin, ray_direction)
    loc_1 = ray_cast_result[1]

    bpy.ops.object.empty_add(location=loc_1)
    bpy.context.object.name = "ray_hit"

    # Local x rotation set to 90 vector
    cam_rot_corrected = Euler(cam_obj.rotation_euler)
    cam_rot_corrected.x = math.radians(90.0)
    cam_rot_corrected = cam_rot_corrected.to_matrix()

    ray_direction = cam_rot_corrected @ Vector((0, 0, -1))

    # Cast the ray from the camera origin onto the object
    ray_cast_result = view_plane.ray_cast(ray_origin, ray_direction)

    ray_cast_result = view_plane.ray_cast(ray_origin, ray_direction)
    loc_2 = ray_cast_result[1]

    bpy.ops.object.empty_add(location=loc_2)
    bpy.context.object.name = "ray_hit"

    # ==Calculate y-shift==
    shift_length = abs(loc_2.z - loc_1.z)
    shift_sign = -1.0 if loc_2.z > loc_1.z else 1.0

    # NOTE: Since the sensor is square and the aspect ratio is 16:9 (1920x1080), the width of the
    #       view frame should correspond to the full length of the sensor.
    # TODO: I think this is probably the broken area.
    #       the sensor and the view frame.
    sensor_length = (frame_points[1] - frame_points[3]).length

    shift_amount = shift_sign * (shift_length / sensor_length)

    # ==Apply y-shift and remove cam local x rotation
    cam_obj.rotation_euler.x = math.radians(90.0)
    cam_data.shift_y = shift_amount

    # ==Cleanup==
    bpy.ops.object.select_all(action="DESELECT")
    view_plane.select_set(True)
    bpy.ops.object.delete(
        use_global=False,
        confirm=False,
    )
```

Answer 1

This can be solved by using the following two things:

• There is a function that converts 3D coordinate values to coordinates on the camera's image plane. The camera coordinates are normalized by the width and height of the plane. bpy_extras.object_utils.world_to_camera_view
• The value of a shift lens is normalized by either the width or height of the image plane. Which one is used depends on the value of cam.data.sensor_fit and the aspect ratio. See my code for more details.

import bpy
import bpy_extras

context = bpy.context
scene  = context.scene
render = scene.render


def auto_remove_camera_vertical_perspective(cam, target):
    ''' target: 3d world coordinates to align with the center of camera coordinates '''

    # Detect deviation of target position from camera center and use it as the shift value
    cam.data.shift_y = 0
    target_ndc = bpy_extras.object_utils.world_to_camera_view(
                                    scene, 
                                    cam, 
                                    target
                                    )
    cam.data.shift_y = target_ndc.y - 0.5

    # Correct if the shift value is normalized by the sensor width of the camera
    aspect_ratio = (render.resolution_x / render.resolution_y 
                    * render.pixel_aspect_x / render.pixel_aspect_y)
    if cam.data.sensor_fit == 'HORIZONTAL' or (
                cam.data.sensor_fit == 'AUTO' and aspect_ratio > 1):
        cam.data.shift_y /= aspect_ratio




cam = scene.camera
target = scene.objects['Empty']

auto_remove_camera_vertical_perspective(cam, target.matrix_world.translation)

Answer 2

Using elements from @tetii's suggestion, this is what I ended up with that appears to work for my use case:

def auto_remove_camera_vertical_perspective_2(cam_obj: bpy_types.Object) -> None:
    context = bpy.context
    scene = context.scene

    # ==Ray cast from camera to get current point in center of camera view==
    # Get the camera's location and rotation
    ray_origin = cam_obj.location
    cam_rot = cam_obj.rotation_euler.to_matrix()

    # Calculate the ray direction from camera's facing
    ray_direction = cam_rot @ Vector((0, 0, -1))

    # Cast the ray from the camera origin into the scene
    depsgraph = bpy.context.evaluated_depsgraph_get()
    ray_cast_result = scene.ray_cast(depsgraph, ray_origin, ray_direction)
    ray_hit_location = ray_cast_result[1]

    # ==Modify camera to remove vertical perspective==
    # Set local x rotation to be parallel to the ground plane
    cam_obj.rotation_euler.x = math.radians(90.0)
    depsgraph.update()

    # Calculate shift y
    cam_space_coords = bpy_extras.object_utils.world_to_camera_view(
        scene,
        cam_obj,
        ray_hit_location,
    )
    shift_y = cam_space_coords.y - 0.5

    # Correct if the shift value is normalized by the sensor width of the camera
    aspect_ratio = (
        scene.render.resolution_x
        / scene.render.resolution_y
        * scene.render.pixel_aspect_x
        / scene.render.pixel_aspect_y
    )

    if cam_obj.data.sensor_fit == "HORIZONTAL" or (
        cam_obj.data.sensor_fit == "AUTO" and aspect_ratio > 1
    ):
        shift_y /= aspect_ratio

    cam_obj.data.shift_y = shift_y
    depsgraph.update()
```