How can we mechanically rotate a branch-cut of a complex function? I'm looking for an algorithmic and mechanical way to define the branch cuts of an arbitrary (or not-so arbitrary, depending on how doable this is) complex function.  Given a function $f(z)$ that has a branch cut along, say, $[1,\infty)$ along the real line, is there a way to algebraically generate a function whose branch cut runs along $1+ci$ (for $c>0$), for instance?  The purpose of this question is to be able to plot arbitrary complex functions with arbitrary choices of branch cut in a computer algebra system like Mathematica, but I believe that this is a math question.  Implementing the answer in Mathematica would be straight-forward once I know how to do this.
As a simple example, can we algebraically transform the function
$$
\sqrt{z^2+1} = \sqrt{|z^2+1|}\exp\left(\frac{1}{2} \operatorname{Arg}({z^2+1}) \right),
$$
which has branch cuts that are rays from $\pm i$ going vertically to infinity, to one that has branch cuts that are rays from $\pm i$ going horizontally to infinity.  I think if I can figure this one out, I'll have what I need.
By way of expaning on my ideas, let's consider the following two situations.
Branch cuts of $f(z) = \sqrt{z^2+1}$ via a Möbius transformation
Two standard ways to define the branch cut of this function are to take the finite segment $[-i,i]$ (defined as the straight-line segment joining $-i$ and $i$) or to take the rays $[i,i\infty)$ and $(-i\infty,-i]$ (defined as the rays going vertically). If we define this function via the complex logarithm as
$$
f(z) = \exp\left(\frac{1}{2}\log\left(z^2+1\right)\right),
$$
and take the standard choice of brach cut for the logarithm, i.e.,
$$
\log(z) = \ln(|z|)+i \operatorname{Arg}(z),
$$
where $-\pi\leq \operatorname{Arg}(z) < \pi$, then this corresponds to the latter choice of  branch cuts ($[i,\infty)$ and $(-\infty,-i]$).
As explained in this answer, we can use a Möbius transformation
$$
z\mapsto \frac{z+i}{z-i},
$$
under which the branch cuts are essentially transformed into each other, and define the function
$$
g(z)=(z-i)\exp\left(\frac{1}{2}\log\left(\frac{z+i}{z-i}\right)\right),
$$
which matches $f(z)$ for $\operatorname{Re}(z)>0$.  We can verify this with Mathematica by plotting these functions (below).
As far as I've been able to work out, this is the only Möbius transformation that works for this case, i.e., it's the only such transformation that makes $f$ and $g$ both be "square roots" of $z^2+1$.

Branch cuts of $\sqrt{z}$ implemented in Mathematica
In Mathematica, $sqrt(z)$ is implemented with a branch cut along the negative real axis, consistent with the choice of branch cut for the logarithm mentioned above.  Indeed, defining a function
f[z_] = Sqrt[Abs[z]] Exp[I*Arg[z]/2]

yields the same results as Sqrt[z], since Mathematica's Arg function restricts its output to be between $-\pi$ and $\pi$.  It is straight-forward to write our own function that rotates the branch cut, as seen here:
sqrt[z_, s_] := Sqrt[Abs[z]] Exp[I/2 Arg[z Exp[-I (s + π)]] + s + π]

This subtracts a phase from the argument $z$, evaluates the Arg function, and then adds the phase back in. Here's an example where we rotate the branch cut by $\pi/4$ (plotting the real part of the function):

Attempt at transforming $\sqrt{z^2+1}$
The same trick doesn't work here.  First of all, we can define this function via the square root function with the different branch cut in the following way.  Define the square root with a rotated branch cut rotated by $\sigma$ as
$$
\operatorname{sqrt}(z) = \sqrt{|z|}\exp\left(
\frac{i}{2}\left(\operatorname{Arg}(ze^{-i(\sigma+\pi)}) + \sigma + \pi\right)
\right).
$$
Then, we can just define
$$
g(z) = \operatorname{sqrt}(z^2+1)
$$
in order to change the branch cut.  Unfortunately, this doesn't rotate the branch cut but rather deforms it, as can be seen in the following plot of the real parts of $\sqrt{z^2+1}$ and $g(z)$ with $\sigma=3\pi/4$ (it starts off linear in the vicinity of the branch point, but then curves away):

So instead, I tried to do something similar to what we did with the square root function directly.  Basically, since
$$
\sqrt{z^2+1} = \sqrt{|z^2+1|}\exp\left(\frac{i}{2}\operatorname{Arg}(z^2+1)\right),
$$
it seems like we should be able to redefine the argument directly in such a way that the branch cut becomes linear in some arbitrary direction by doing something like
$$
\operatorname{Arg}(z^2+1) = \operatorname{Arg}\left(\left(ze^{-i(\sigma+\pi)}\right)^2+1\right) + (?),
$$
but I cannot figure out what $\sigma$ should be.
Is there a general recipe for this?
 A: Here is a partial answer that works for the case of $f(z)=\sqrt{z^2+1}$ and will likely work for general power functions.  I don't show that it will work for general functions, but I think so. In addition, I haven't determined generally what needs to happen in order to plot the "entire" function, by which I mean the function over all Riemann sheets, without overlap.  I've figured it out for the function $f$ by trial and error, but I'm not sure about the careful mathematics.
In addition, I don't know how to make this work for functions that aren't expressible as power functions of polynomials like $f$.

The trick is to play fast and loose with the complex arithmetic of root rules and just factor the function into "linear" factors, i.e.,
$$
\sqrt{z^2+1} "=" \sqrt{z+i}\sqrt{z-i}.
$$
These two expressions are not actually equal under the standard choice of branch cut for the square root, as can be seen by plotting their real parts (the two functions are identical everywhere except where you can see the orange):

The branch cuts have been rotated by $\pi/2$ in opposite directions.  In order to rotate the branch cuts arbitrarily, we use the $\operatorname{sqrt}$ function from the OP:
$$
\operatorname{sqrt}(z;\sigma) = \sqrt{|z|}\exp\left(
\frac{i}{2}\left(\operatorname{Arg}(ze^{-i(\sigma+\pi)}) + \sigma + \pi\right)
\right).
$$
This allows us to choose the branch cuts for each square root independently.  Note that taking $\sigma=-\pi$ gets us the standard choice for the branch cut of $\sqrt{\cdot}$, i.e.,
$$
\sqrt{z} = \operatorname{sqrt}(z;-\pi).
$$
By experimenting, we can determine that one choice of branch cuts makes the expressions above "equal", i.e.
$$
\sqrt{z^2+1} = \operatorname{sqrt}\left(z^2+1;-\pi\right) = \operatorname{sqrt}\left(z+i;-\pi-\frac{\pi}{2}\right)
\operatorname{sqrt}\left(z-i;-\pi+\frac{\pi}{2}\right),
$$
which can likely be verified algebraically, but here's a plot:

Now, we can rotate the two branch cuts independently of each other. Below is a plot that includes two Riemann sheets for this function.  The trick to making two sheets here is to keep one of the square roots "on the first sheet" and move the other up one sheet.  Here's what I mean.  We'll plot these two functions:
\begin{align}
g_{\textrm{yellow}}(z)&=
\operatorname{sqrt}\left(z+i;-\pi-\frac{\pi}{4}\right)
\operatorname{sqrt}\left(z-i;-\pi+\frac{\pi}{4}\right)\\
g_{\textrm{blue}}(z)&=
\operatorname{sqrt}\left(z+i;\pi-\frac{\pi}{4}\right)
\operatorname{sqrt}\left(z-i;-\pi+\frac{\pi}{4}\right).
\end{align}
Notice that the only difference between the two is that the argument of one of the $\operatorname{sqrt}$ functions has increased by $2\pi$ to $\pi-\frac{\pi}{4}$.  Here is the result:

