So here's the better way:
new_x = old_x * cos(angle) - old_y * sin(angle) new_y = old_x * sin(angle) + old_y * cos(angle)
With that out of the way, what your code is doing is composing three rotations around axes. You haven't explained why you think the results are incorrect. Chances are your intuitions about rotations are a bit off and the code is working just fine.
So what's an example of inputs into this function that produce output that you consider incorrect? What result did you expect instead?