Let $q$ the position of your robotic arm motors, so its value describe completely the position of your manipulator at each time. First we suppose that the motors of our robot have no limits so, abusing of the notation, we write $q \in \mathcal Q$, i.e., the configuration space and the joint's velocities $ \dot q \in \mathbb R^n$ for all $\mathcal Q$. Now we call $p = [\matrix{x^\top \psi^\top}]^\top$ the column vector with position and the euler's angles. Now, the differential direct kinematics reads

$$ \dot p = J(q) \;\dot q $$

Where the jacobian $J$ has to parts $$ J=\pmatrix{J_x \cr J_\psi} $$

This equation can be readed as "for infinitesimal changes in the joint's positions, give me the infinitesimal changes in the end-effector's positions/orientation". In fact velocity deals with infinitesimal displacement in infinitesimal time.

The limits of the operations space $\Omega \in \mathbb R^3$ of a manipulator are the points where no infinitesimal changes of the joints positions gives infinitesimal changes in space positions. In other words, at these points it is not possible to achieve displacements in, at least, one direction. This means that it the joint configurations $q_0$ is in the limit of hte workspace if

$$

$$

Each motor have its limit, so for each ith motor $q_i \in [q_{i}^0, q_i^1]$ , so for simplicity we write that $q \in \mathcal Q \subset \mathbb R^n$ for an $n$ dof robot.

I will explain how to determinate the "normal" three dimensional workspace and then extrapolate to orientation. I hope to be useful.

Let $q$ the position of your robotic arm motors, so its value describe completely the position of your manipulator at each time. First we suppose that the motors of our robot have no limits so, abusing of the notation, we write $q \in \mathcal Q$, i.e., the configuration space and the joint's velocities $ \dot q \in \mathbb R^n$ for all $\mathcal Q$. Now we call $p = [\matrix{x^\top \psi^\top}]^\top$ the column vector with position and the euler's angles. Now, the differential direct kinematics reads

$$ \dot p = J(q) \;\dot q $$

Where the jacobian $J$ has to parts $$ J=\pmatrix{J_x \cr J_\psi} $$

This equation can be readed as "for infinitesimal changes in the joint's positions, give me the infinitesimal changes in the end-effector's positions/orientation". In fact velocity deals with infinitesimal displacement in infinitesimal time.

The limits of the operations space $\Omega \in \mathbb R^3$ of a manipulator are the points where no infinitesimal changes of the joints positions gives infinitesimal changes in space positions. In other words, at these points it is not possible to achieve displacements in, at least, one direction. This means that it the joint configurations $\bar{q}$ is in the limit of hte workspace if

$$ \mbox{Rank}\;( J_x(\bar{q})) < 3 $$

because there is at least one direction of the three dimensional space $\mathbb R^3$ that is outside the co-domain of the linear mapping given by the matrix $J_x(\bar q)$

The same reasoning can be applied to orientation. We looks for points in $\mathcal Q$ where infinitesimal changes in the joint's positions give no change in some orientation parameter, say an euler angle. The problem of using the previous formula as $\mbox{Rank}\; J_\psi$ is that the matrix $J_\psi$ could have a rank less than 3 because a representations singularity. That can be overrided using another minimal orientation representation like quaternions, or simply changing the used euler angles in the singular points.

Now if each motor have its limit we write $q_i \in [q_{i}^0, q_i^1]$. We call $\partial \mathcal Q $ the set of all configurations $q$ for which at least one component is in its limit. We have that for $ q \in \partial \mathcal Q $ there is at least one $i$ for which $\dot q_i \in \mathbb R^+$ XOR $\dot q_i \in \mathbb R^-$, i.e., for the joint $q_i$, when it is at its upper limit, necessarily $\dot q \le 0$, and, when it is at its lower limit, necessarily $\dot q \ge 0$. For the limits of the workspace $\Omega \subset \mathbb R^3$ we have that the formula $\mbox{Rank}\;J_x < 3$ cannot be used in $\partial \Omega$, because at these points the domain of $J_x$ is not a vector space. At the points $\bar q \in \partial \mathcal Q$ there is one direction that is forbidden for $\dot {\bar q}$, the mapping of that directions $J_x \dot{\bar q}$ is also a forbidden direction if and only if it is linearly independent of the mapping of all possible velocities $J_x \dot q_j$ in the three-dimensional space, then it is a limit of the workspace. In other words, a forbidden joint velocity maps into a forbidden end-effector displacement if and only if, there is not another allowed infinitesimal joint displacement that maps into a displacement in that direction. The same reasoning can be applied to orientation, taking into account the representations singularities.

So, your software shout do this:

1. in $\mathcal Q$ just seek for the values of $q$ where $\mbox{Rank}\;J_\psi < 3$ taking into account the representation singularities.
2. for each point in $\partial \mathcal Q$ you have to identify the forbidden velocity, say $\dot z$ and take its mapping $v = J_\psi \dot z$. Now construct the matrix $J'_\psi$ removing the $i$ column of each joint that lies in its limits. If $v^\top (J'_\psi) =0 $ then the direction is a limit in the orientation. Always taking into account the representations singularities!