The Extract Transform VOP can extract either one of translate, rotate or scale from a 4x4 matrix. This is how I extract the rotate out of the transform matrix in the example above.
If your question is given two different sets of the same geometry in different orientations and either has the same shape and topology, extract the rotate values. Well that isn't so straightforward but quite possible to do but since you don't have the advantage of the transform matrix that did the transform, you will run in to issues, especially if you wish to calculate the rotate with animated targets that roll more than 360 degrees (or 180). Solutions are prone to flipping about the 360 (and 180) degree values and this is to be expected. Whenever the two vectors are parallel regardless of where they start and end.
Having said all of that, it's actually pretty easy. For the cases where both sets of geometry are not deforming, you can take any point from the two sets of geometry, calculate their surface N normals, then in VOPs use the Align VOP (using the dihedral() function, wikipedia it or just accept it for what it does) and this generates a 3x3 matrix that holds the rotation to align the first vector to the second vector. Then use the Extract Transform VOP again to extract three Euler rotate values. The sweet part is that VOPs automatically type casts the 3x3 matrix to a 4x4 matrix correctly for you without having to insert a VOP to do it. Part of the Ease-Of-Use drive over the years.
As for the flipping case, what happens when the two vectors are perfectly parallel (dihedral returns 0) compared to when they are almost parallel on either side of perfect alignment. One will go from 0.0000001 to 0. The other from 359.9999 to 0. The latter gives you the flip and as you approach 0, the significant digits will easily exceed the precision available and that goes for doubles as well. Actually as you approach 0 from the 359 side of the equation, you reach infinity (indeterminate). Don't worry. I'm not asking for proofs. Just need to be aware of this. This is the same issue with look-ats and character kinematics and joint flipping. Yes you can use the Rotate VOP (dihedral()) to calculate look-at rotations.
You could also construct a quaternion in a 4x4 matrix to do the rotates but even then it will still be prone to flipping around the point where your two vectors sweep through being parallel.
What do I do to remove this ambiguity? Do everything in Houdini why of course and keep everything procedural.
Another way is to use CHOPs and works well for the case where the geometry is animated and deforming. Record your singe set of deforming animated geometry with a Geometry CHOP. You will calculate the transforms relative to the starting frame geometry position. This also includes transforms and rotates. You will see discontinuities in the rotate channels. Just append a Transform CHOP and this will properly accumulate the rotates and remove the discontinuities in the rotate channels. Export the rotates to wherever you want. At this point it's a toss as to what's faster: point() exp
ression to extract rotates in to paris or CHOPs...
CHOPs also has ways of parsing vectors (t channels, r channels) and can do full matrix transforms if the t, r and s channels are there for every point. As usual, the matrices are handled for you so you just need to know that 4x4 matrices can hold Transforms. No magic what so ever.
That is what I would do if it was mission critical and live with the hit you take with CHOPs.