import bpy
import math
import numpy
import mathutils

# @brief オブジェクトを指定座標へ向けるための回転軸と回転角を計算する
# @param [in] オブジェクト
# @param [in] 向ける方向（座標）
# @return 回転軸,回転角(ラジアン)
def calc_lookat_param(obj,lookat):
    
    # オブジェクトのローカルのZ軸のベクトルを取得
    mat = obj.matrix_world
    localZ = mathutils.Vector((mat[0][2],mat[1][2],mat[2][2]))
    
    va = mathutils.Vector(lookat) - obj.location
    vb = mathutils.Vector(localZ)
    
    va.normalize()
    vb.normalize()
    
    # 外積
    axis = mathutils.Vector.cross(va,vb)
    axis.normalize()
    
    # 内積
    th = mathutils.Vector.dot(va,vb)
    
    # 角度算出
    rad = -math.acos( numpy.clip(th,-1.0,1.0) )

    return axis , rad
    

    
# @brief 任意軸回転
# @param [in,out] obj 回転するオブジェクト
# @param [in] axis 回転軸
# @param [in] radian 回転角をラジアンで指定
# @return なし
# @sa https://www.study.suzulang.com/bpy-calculation/bpy-arbitrary-axis-rotation
def rotate_object(obj,axis,radian):
    
    
    rot_mat= mathutils.Matrix.Rotation( radian, 4, axis ) 
    
    # decompose world_matrix's components, and from them assemble 4x4 matrices
    orig_loc, orig_rot, orig_scale = obj.matrix_world.decompose()
    #
    orig_loc_mat   = mathutils.Matrix.Translation(orig_loc)
    orig_rot_mat   = orig_rot.to_matrix().to_4x4()
    orig_scale_mat = (mathutils.Matrix.Scale(orig_scale[0],4,(1,0,0)) @ 
                      mathutils.Matrix.Scale(orig_scale[1],4,(0,1,0)) @ 
                      mathutils.Matrix.Scale(orig_scale[2],4,(0,0,1)))
    #
    # assemble the new matrix
    obj.matrix_world = orig_loc_mat @ rot_mat @ orig_rot_mat @ orig_scale_mat 
    
    
        
    

# 使用例
obj = bpy.context.active_object

axis , rad = calc_lookat_param( 
    obj,
    bpy.context.scene.cursor.location)

# 実際に回転
rotate_object(obj,axis,rad)