Hi,
If the XML file contains the MeshValueCollection in the 'right' format, then you can load in as a mesh function in python:
  mesh=Mesh("mesh_file_name")
  mvc = MeshValueCollection("size_t", mesh, "test_mesh_bcs.xml", 1)
Then you need convert to mesh function.  The problem is how to find the global index of facets based on the local index. 
  mvc.values()   # this would return a dictionary to indicate the local index and values
I know there is function in DOLFIN to find the connection between local index and global index for vertices in the GENERICDOFMAP class, but i do not know if there exists such a function to find the map between local index and global one for facets. 
But I think you can find that through the following way: (mvs.values() contains the cell index and local index in that cell)
- Based on the cell index,  find the connectivity of the cell with the facets. Then with that connectivity, we get the global index of facets for each local one 
- Create the mesh function on facets and then we just need to assignment the corresponding value to finish the construction of the mesh function you want. 
Hope it is not a mess. And if anyone better way to find the connection between local and global index for facets, that would be great and be beneficial to me also.