As some people on tex.sx, I used the tikz-3dplot and changed some commands to use Tait-Bryan convention instead of the default Euler convention.
As you can see on the following picture, the rotation works fine and the rotated xy-plane (magenta one) and the xz-plane (the black one) are drawn properly (each of them, centered in (0,0,0) touches the top of their two respective unit vectors.
Nevertheless, the third plane, yz (the red one) is not on the good plane.
The workaround to draw this plane (the cyan one) was to make a 3d scope.
%Using directly the 3d library instead
\tdplotsetrotatedcoords{\yaw}{\pitch}{\roll}
\begin{scope}[tdplot_rotated_coords,canvas is yz plane at x=0]
\draw[dashed,cyan,->,thick] (0:1) arc (0:350:1);
\end{scope}
This clearly shows that the rotation does its work, it's in the yz plane and the 0 stands on the y axis. The rotation is in the good direction.
So, this might be because of the redefinition of the thetaplane command :
%Instead of permuting the coordinates as does Jeff, I rotate the planes as in the previous command.
\renewcommand{\tdplotsetrotatedthetaplanecoords}[1]{%
\tdplotresetrotatedcoordsorigin
\tdplotsetrotatedcoords{\tdplotalpha + #1}{\tdplotbeta}{\tdplotgamma+90}
}
So, am I doing something wrong ?
Anyway, there is my code with all the renewcommand :
\documentclass{standalone}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsthm}
\usepackage[usenames,dvipsnames]{xcolor}
\usepackage{tikz}
\usepackage{tikz-3dplot}
\usepackage{pgfplots}
\usepgfplotslibrary{external}
%You get back to Euler convention by commenting out the following lines.
%Redefine all the tikz-3dplot according to Tait-Bryan convention
%This performs the calculation to define the main coordinate frame orientation style, and is also used to transform a coordinate from the main coordinate frame the the screen coordinate frame
%In comparison to genuine tikz-3dplot, this changes the orientation, so z points downwards when turnen around x-axis by 90 degres
\renewcommand{\tdplotcalctransformmainscreen}{%
%
\tdplotsinandcos{\sintheta}{\costheta}{\tdplotmaintheta}%
\tdplotsinandcos{\sinphi}{\cosphi}{\tdplotmainphi}%
%
\tdplotmult{\stsp}{\sintheta}{\sinphi}%
\tdplotmult{\stcp}{\sintheta}{\cosphi}%
\tdplotmult{\ctsp}{\costheta}{\sinphi}%
\tdplotmult{\ctcp}{\costheta}{\cosphi}%
%
%determine rotation matrix elements for display transformation
%it's line major meaning that rab_ is the first line, second column
\pgfmathsetmacro{\raarot}{\cosphi}%
\pgfmathsetmacro{\rabrot}{-\sinphi}%
%NOTE: \rac is zero for this rotation, where z^c always points vertical on the page
\pgfmathsetmacro{\racrot}{0}%
\pgfmathsetmacro{\rbarot}{\ctsp}%
\pgfmathsetmacro{\rbbrot}{\ctcp}%
\pgfmathsetmacro{\rbcrot}{-\sintheta}%
%NOTE: third row of rotation matrix not needed for display since screen z is always flat on the page. It is, however, needed when performing further transformations using the Euler transformation.
\pgfmathsetmacro{\rcarot}{\stsp}%
\pgfmathsetmacro{\rcbrot}{\stcp}%
\pgfmathsetmacro{\rccrot}{\costheta}%
%
}
%determines the rotation matrix for transformation from the rotation coordinate frame to the main coordinate frame. This also defines the rotation to produce the rotated coordinate frame.
%In comparison to genuine tikz-3dplot, it does a zyx rotation sequence according to Tait-Bryan convention.
\renewcommand{\tdplotcalctransformrotmain}{%
%perform some trig for the Euler transformation
\tdplotsinandcos{\sinalpha}{\cosalpha}{\tdplotalpha}
\tdplotsinandcos{\sinbeta}{\cosbeta}{\tdplotbeta}
\tdplotsinandcos{\singamma}{\cosgamma}{\tdplotgamma}
%
\tdplotmult{\sasb}{\sinalpha}{\sinbeta}
\tdplotmult{\sasg}{\sinalpha}{\singamma}
\tdplotmult{\sasbsg}{\sasb}{\singamma}
%
\tdplotmult{\sacb}{\sinalpha}{\cosbeta}
\tdplotmult{\sacg}{\sinalpha}{\cosgamma}
\tdplotmult{\sasbcg}{\sasb}{\cosgamma}
%
\tdplotmult{\casb}{\cosalpha}{\sinbeta}
\tdplotmult{\cacb}{\cosalpha}{\cosbeta}
\tdplotmult{\cacg}{\cosalpha}{\cosgamma}
\tdplotmult{\casg}{\cosalpha}{\singamma}
%
\tdplotmult{\cbsg}{\cosbeta}{\singamma}
\tdplotmult{\cbcg}{\cosbeta}{\cosgamma}
%
\tdplotmult{\casbsg}{\casb}{\singamma}
\tdplotmult{\casbcg}{\casb}{\cosgamma}
%
%determine rotation matrix elements for Euler transformation
\pgfmathsetmacro{\raaeul}{\cacb}
\pgfmathsetmacro{\rabeul}{\casbsg - \sacg}
\pgfmathsetmacro{\raceul}{\sasg + \casbcg}
\pgfmathsetmacro{\rbaeul}{\sacb}
\pgfmathsetmacro{\rbbeul}{\sasbsg + \cacg}
\pgfmathsetmacro{\rbceul}{\sasbcg - \casg}
\pgfmathsetmacro{\rcaeul}{-\sinbeta}
\pgfmathsetmacro{\rcbeul}{\cbsg}
\pgfmathsetmacro{\rcceul}{\cbcg}
%DEBUG: display euler matrix elements
%\raaeul\ \rabeul\ \raceul
%
%\rbaeul\ \rbbeul\ \rbceul
%
%\rcaeul\ \rcbeul\ \rcceul
}
%determines the rotation matrix for transformation from the main coordinate frame to the rotated coordinate frame.
%In comparison to genuine tikz-3dplot, it does a zyx rotation sequence according to Tait-Bryan convention.
%It's the transposed of the matrix calculated in tdplotcalctransformrotmain
\renewcommand{\tdplotcalctransformmainrot}{%
%perform some trig for the Euler transformation
\tdplotsinandcos{\sinalpha}{\cosalpha}{\tdplotalpha}
\tdplotsinandcos{\sinbeta}{\cosbeta}{\tdplotbeta}
\tdplotsinandcos{\singamma}{\cosgamma}{\tdplotgamma}
%
\tdplotmult{\sasb}{\sinalpha}{\sinbeta}
\tdplotmult{\sasg}{\sinalpha}{\singamma}
\tdplotmult{\sasbsg}{\sasb}{\singamma}
%
\tdplotmult{\sacb}{\sinalpha}{\cosbeta}
\tdplotmult{\sacg}{\sinalpha}{\cosgamma}
\tdplotmult{\sasbcg}{\sasb}{\cosgamma}
%
\tdplotmult{\casb}{\cosalpha}{\sinbeta}
\tdplotmult{\cacb}{\cosalpha}{\cosbeta}
\tdplotmult{\cacg}{\cosalpha}{\cosgamma}
\tdplotmult{\casg}{\cosalpha}{\singamma}
%
\tdplotmult{\cbsg}{\cosbeta}{\singamma}
\tdplotmult{\cbcg}{\cosbeta}{\cosgamma}
%
\tdplotmult{\casbsg}{\casb}{\singamma}
\tdplotmult{\casbcg}{\casb}{\cosgamma}
%
%determine rotation matrix elements for Euler transformation
\pgfmathsetmacro{\raaeul}{\cacb}
\pgfmathsetmacro{\rabeul}{\sacb}
\pgfmathsetmacro{\raceul}{-\sinbeta}
\pgfmathsetmacro{\rbaeul}{\casbsg - \sacg}
\pgfmathsetmacro{\rbbeul}{\sasbsg + \cacg}
\pgfmathsetmacro{\rbceul}{\cbsg}
\pgfmathsetmacro{\rcaeul}{\sasg + \casbcg}
\pgfmathsetmacro{\rcbeul}{\sasbcg - \casg}
\pgfmathsetmacro{\rcceul}{\cbcg}
%
%DEBUG: display euler matrix elements
%\raaeul\ \rabeul\ \raceul
%
%\rbaeul\ \rbbeul\ \rbceul
%
%\rcaeul\ \rcbeul\ \rcceul
}
%\tdplotsetthetaplanecoords{\phi}
%this places the rotated coordinate system such that it's x'-y' plane coincides with a "theta plane" for the main coordinate system: This plane contains the z axis, and lies at angle \phi from the x axis.
%#1: user-specified \phi angle from x-axis
%Since the rotation sequence changed, we also have to redefine the theta plane. We basically turns the zx-plane to the theta plane and flip the xy plane of 90° to make him coincide with the theta plane.
\renewcommand{\tdplotsetthetaplanecoords}[1]{%
%
\tdplotresetrotatedcoordsorigin
\tdplotsetrotatedcoords{#1}{0}{90}%
}
%%\tdplotsetrotatedthetaplanecoords{\phi'}
%%this places the rotated coordinate system into the "theta plane" for the current rotated coordinate system, at user-specified angle \phi'. Note that it replaces the current rotated coordinate system
%%#1: user-specified \phi' angle from x'-axis
%Instead of permuting the coordinates as does Jeff, I rotate the planes as in the previous command.
\renewcommand{\tdplotsetrotatedthetaplanecoords}[1]{%
\tdplotresetrotatedcoordsorigin
\tdplotsetrotatedcoords{\tdplotalpha + #1}{\tdplotbeta}{\tdplotgamma+90}%
% \message{DEBUG theta plane coords : #1 : \racrc, \rbcrc}
}
\begin{document}
\def\roll{30}
\def\pitch{50}
\def\yaw{40}
\def\xMainRot{100}
\def\zMainRot{30}
\tdplotsetmaincoords{\xMainRot}{\zMainRot}
\begin{tikzpicture}[scale=4,tdplot_main_coords,every node/.append style={transform shape}]
%
\draw[thick,->] (0,0,0) -- (1,0,0) node[anchor=north east]{$x$};
\draw[thick,->] (0,0,0) -- (0,1,0) node[anchor=north west]{$y$};
\draw[thick,->] (0,0,0) -- (0,0,1) node[anchor=south east]{$z$};
%
\tdplotsetrotatedcoords{\yaw}{\pitch}{\roll}
\draw[thick,color=magenta,tdplot_rotated_coords,->] (0,0,0) -- (1,0,0) node[anchor=north]{$x$};
\draw[thick,color=magenta,tdplot_rotated_coords,->] (0,0,0) -- (0,1,0) node[anchor=west]{$y$};
\draw[thick,color=magenta,tdplot_rotated_coords,->] (0,0,0) -- (0,0,1) node[anchor=south west]{$z$};
%
%Drawing the three planes of the new frame
\tdplotdrawarc[tdplot_rotated_coords,dashed,color=magenta]{(0,0,0)}{1}{0}{360}{anchor=north,transform shape}{}
\tdplotsetrotatedthetaplanecoords{0}
\tdplotdrawarc[tdplot_rotated_coords,dashed,color=black]{(0,0,0)}{1}{0}{360}{transform shape}{}
\tdplotsetrotatedcoords{\yaw}{\pitch}{\roll}
\tdplotsetrotatedthetaplanecoords{90}
\tdplotdrawarc[tdplot_rotated_coords,dashed,color=red]{(0,0,0)}{1}{0}{360}{transform shape}{}
%
%Using directly the 3d library instead
\tdplotsetrotatedcoords{\yaw}{\pitch}{\roll}
\begin{scope}[tdplot_rotated_coords,canvas is yz plane at x=0]
\draw[dashed,cyan,->,thick] (0:1) arc (0:350:1);
\end{scope}
%NOTE: the rotation does its work, it's in the yz plane and the 0 stands on the y axis and the rotation is in the good direction.
\end{tikzpicture}
\end{document}