Could be related to: pull back on quadrilateral not working if gdim = 3 · Issue #2109 · FEniCS/dolfinx · GitHub
which I will have a look at later today.
Could you print mesh.geometry.dim. I bet it is 3. If you prune the z dimension from your points in
i.e.
points=mesh.points[:,:2]
I think it should work.