How to get the updated rotation of an object during soft body physics simulation?
Issue with all vert coords in vg summing to zero.
Made a quick demo to have a squiz,
Using 2.91 ran into this issue
Default Cube set up to drop on plane using physics. An all vertex group &quot;Group&quot; is added, and two constraints, as demonstrated with code above.
All constraints use default values unless otherwise stated.
Obvservations:
For copy rotation constraint.
Using no subtarget --&gt; Empty same rotation as cube.
When using all vert vg &quot;Group&quot; as a subtarget  --&gt; Resulted in Identity (zero) rotation.
Using a group with only 3 vertices --&gt; Empty rotates, but is aligned to normal of 3 verts plane.
Editing the cube by moving one corner such that it is no longer symmetrical also sees rotation when all verts used as subtarget.
Similarly if both constraints are replaced with one copy transform constraint, in which case &quot;Group&quot;  makes the empty disappear (scale to zero). Despite there being only unit scale on all objects involved
Didn't try but the other similar constraint would be parent to.
Conclusion.
Either don't use a subtarget, or make sure the sum of all vert coords  in  subtarget group do not sum to zero,  for constraints that involve rotation.
Really  only need the subtarget if also deforming the mesh.
Perhaps someone can elaborate whether this behaviour re using all verts vg as a subgroup not rotating is by design, or a bug of sorts, when caused by the sum of all cube verts being zero a zero vector...
With one of the working options selected, used script below to add an empty at each tenth frame with same matrix world as constrained empty.
import bpy
scene = bpy.context.scene
fs, fe = scene.frame_start, scene.frame_end    
tg = bpy.data.objects[&quot;Empty&quot;]

for frame in range(fs, fe):
   scene.frame_set(frame)
   print(frame, tg.matrix_world.to_euler('XYZ')) # This gives me completely random rotation vectors
   if not frame % 10:
       bpy.ops.object.empty_add()
       mt = bpy.context.object
       mt.matrix_world = tg.matrix_world
       mt.name = f&quot;frame_{frame}&quot;
       mt.show_name = True

With the physics baked, ran it and the added empties align with the constrained empty at their respective frame.
PS can make an all vert group with
group.add(range(len(mesh.vertices)), 1.0, &quot;ADD&quot;)