How does blender calculate poseBone.matrix? How does blender calculate a poseBone.matrix through the hierarchy chain? 
A bone's final matrix is a transformation matrix which defines the bone in armature-space in the current pose (M1). I want to calculate this from other matrices, given a really simple armature with two bones, the second bone parented to the first.

armature = bpy.data.objects["Armature"]
M0 = armature.pose.bones[0].matrix       #final matrix for the first bone
M1 = armature.pose.bones[1].matrix       #final matrix for the second bone


A bone's bind-pose or rest-pose matrix is its transformation from the armature origin:

l0 = armature.data.bones[0].matrix_local


I can calculate the local rest matrix relative to the bone's parent:

l1 = armature.data.bones[1].matrix_local * l0.copy().inverted()


I can seem to calculate the final matrix of the first bone (M0) like this:

b0 = armature.pose.bones[0].matrix_basis    #local bone transform
l0*b0    #this equals M0


However, if I want to calculate the second bone's matrix, which is parented to the first, I can not get the same M1 matrix back while trying this:

b1 = armature.pose.bones[1].matrix_basis    #local bone transform
l0*b0*l1*b1    #I wanted to go through the parent chain like this but couldn't get back M1


However, the resulting matrix has the rotation part of M1, but the translation part is different.

The reason I want to do this to see how to do a bone hierarchy system for myself, but I seem to be missing something.
Look at this link. the difference between matrices is described here:
How can I manually calculate bpy.types.PoseBone.matrix using Blender&#39;s Python API?