## 2 Methods  
### Algorithm Multiplexing  
The Modified Gerchberg-Saxton algorithm (MGSA) used in this study could transform phase-only function (POF) with optical wavelength or imaging position [12; 13]. The POFs with distinct encryption conditions on the Fresnel transform (FrT) plane were combined as a single POM (phase-only mask). The study applied position multiplexing to realize the multi-depth image CGH display. Using the x-axis as the axis of rotation as an example, the geometric coordinate of Fresnel diffraction between the corresponding imaging tilted plane and input plane [13] is shown in Figure 1. After the image phase information traveled through the Fresnel diffraction, the output image parallel to the reference image would be shown on the output plane. The original field signal placed on the input plane, after Fresnel diffraction, could acquire the diffraction field image \(G(x,y)\) on the output plane, which is parallel to the reference plane. To obtain the distance between a point on the tilted output plane and the corresponding point on the input plane, a calculation, shown as Equation (1), was required.  
\[r_{x}=\sqrt{(x_{0}-x^{\prime})^{2}+(y_{0}-y^{\prime}\cos\theta_{x})^{2}+(z-y^{ \prime}\sin\theta_{x})^{2}} \tag{1}\]  
where \(\theta x\) stands for the angle of \(x\) axial rotation. When the output plane revolves the \(x\) axis and is not parallel to the input plane, only the Fresnel diffraction of the object on the tilted plane [14; 15] would be calculated in order to parallel the diffracted output plane with the reference plane. To acquire the tilted Fresnel equations, the signal \(g(x_{0},y_{0})\) on the input plane is transformed, shown as Equation (2) below:  
\[\begin{split} G_{x}(x^{\prime},y^{\prime})=\exp\bigl{(}i\frac{2 \pi}{\lambda}r_{x}\bigr{)}&\int\limits_{-\infty}^{\infty}g(x_{0},y_{0})\exp\Bigl{[}i\frac{\pi}{\lambda r_{x}}\bigl{(}x_{0}^{2}+y_{0}^{2}\bigr{)} \Bigr{]}\\ &\times\exp\Bigl{\{}-i\frac{2\pi}{\lambda r_{x}}[x_{0}x^{\prime} +y_{0}y^{\prime}\cos\theta_{x}+(z_{0}-r_{x})y^{\prime}\sin\theta_{x}]\Bigr{\}} dx_{0}dy_{0}\end{split} \tag{2}\]  
where \(\lambda\) is the wavelength and \(z_{0}\) is the distance between the center point on the tilted input plane and the center point on the output plane. Assuming \(N\) pieces of the target image being \(g_{n}(x^{\prime},y^{\prime})\), \(n=1\)-\(N\), each image corresponds to the different angles of the inclination \(\theta x\). The individual phase information \(\psi_{\theta x}(x_{0},y_{0})\) was acquired through MGSA calculation.  
When the value above is substituted in Equation (1), it could be simplified as Equation (3), shown below:  
\[\begin{array}{c}G_{x}(x^{\prime},y^{\prime})=\text{FrT}\{\exp\{j\psi_{\theta_{xn }}(x_{0},y_{0})\};\lambda;z_{0};\theta_{xn}\}\\ \qquad\qquad\qquad\qquad=\hat{g}_{\theta_{xn}}(x^{\prime},y^{\prime})\exp\{j \varphi g_{\theta_{xn}}(x^{\prime},y^{\prime})\}\end{array} \tag{3}\]  
where \(\varphi g_{\theta_{xn}}(x^{\prime},y^{\prime})\) is the individual phase signal of \(\psi^{\prime}_{\theta_{xn}}(x_{0},y_{0})\), \(n=1\sim N\) after FrT, shown in Figure 2. Individual phases were integrated into a single piece of phase-only mask (POM), as shown in Equation (4):  
\[\exp\{j\psi^{\theta_{xn}}_{T}(x_{0},y_{0})\}=\exp\{\sum_{n=1}^{N}\exp\{j\psi^{ \prime}_{\theta_{xn}}(x_{0},y_{0})\}\} \tag{4}\]