In here there are functions that allow to plot geodesics in the upper half plane (namely, half circles perpendicular to the x axis and vertical lines). Are there similar functions that do the same thing but only in the fundamental domain of SL2(Z), so that the geodesic is folded up into the domain |x|<12 and x2+y2>1?
What I would like to have is something like the picture in page 5 of this paper.