The main objective of this work is to show that this method can provide the optimal experimental conditions (structure and beam dimensions, angle of incidence) for which the structure operates effectively. On the one hand, the incident beam (Gaussian, Hermite-Gaussian, Laguerre-Gaussian, etc.) is known through its angular spectrum expansion \(\overrightarrow{F}(k_{X},k_{Y},Z)\) given in its frame, namely the \(OXYZ\) frame (see Figure 1) corresponding to an average propagation direction along the \(OZ\) direction. Generally, the angle of incidence on the metasurface is arbitrary. It is named \(\theta_{m}\), as shown in Figure 1.  
On the other hand, the monomode modal method, explained in more detail in [14], allows the determination of the transmission and reflection Jones matrices in the \((TE,TM)\equiv(s,p)\) basis for any illumination direction characterized by the angle of incidence \(\theta\) and the azimuthal angle \(\psi\). For the transmission, we have:  
\[t(\theta,\psi)=\left( \begin{array}{cc}t_{ss}&t_{sp}\\ t_{ps}&t_{pp}\end{array} \right) \tag{1}\]  
This matrix provides the transmission (amplitude and phase) for any incident plane wave illuminating the metasurface. To determine the transmitted beam, we combine the angular spectrum expansion with the Jones formalism through a variable change allowing the expression of all variables in the metasurface frame (\(Oxyz\)) in which the Euler angles (\(\theta\) and \(\psi\)) are defined (see Figure 2). This implies a variable change to express the angular spectrum in this frame. First, we need to express the Jones matrix of each plane wave as a function of its wave-vector components \(\overrightarrow{k}(k_{x},k_{y},k_{z})\) through the relations between the latter and the angles \(\theta\) and \(\psi\) given by:  
\[ \begin{array}{rcl}\theta&=&\cos^{-1}(\frac{k_{z}}{\sqrt{\varepsilon_{1}k_{0} }})\\ \psi&=&\cos^{-1}(\frac{k_{x}}{k_{\parallel}})\end{array}  \tag{2}\]  
with \(k_{z}=\sqrt{\varepsilon_{1}k_{0}^{2}-k_{\parallel}^{2}}\)\(k_{\parallel}\) is the tangential component of the wave vector (\(k_{\parallel}=\sqrt{k_{x}^{2}+k_{y}^{2}}\)) and \(k_{0}=\frac{\omega}{\varepsilon}\) is the wave-vector modulus in vacuo. Figure 1: Schematic of the model with the different frames involved in the theoretical model: \(Oxyz\) the metasurface frame, \(OXYZ\) the incident beam frame, \(OX^{\prime\prime}Y^{\prime\prime}Z^{\prime\prime}\) the reflected beam frame and \(OX^{\prime}Y^{\prime}Z^{\prime}\) the transmitted beam frame. \(\theta_{m}\) is the angle of incidence. Second, a basis change from \((x,y)\rightarrow(s,p)\rightarrow(x,y)\) is necessary to determine the transmitted angular spectrum \(\overrightarrow{F_{t}}(k_{x},k_{y})\) components in the \(Oxyz\) frame of the metasurface:  
\[\overrightarrow{F_{t}}(k_{x},k_{y})=\wp^{t}\cdot t(\theta,\psi)\cdot\wp\cdot \overrightarrow{F}(k_{x},k_{y}) \tag{3}\]  
The \({}^{t}\) operator is the transpose matrix operator, and \(\wp\) is the basis change matrix \((s,p)\rightarrow(x,y,z)\) given by:  
\[\wp(k_{x},k_{y})=\left( \begin{array}{cc}-\dfrac{k_{x}k_{z}}{ \sqrt{\varepsilon_{1}}k_{0}k_{\parallel}}&\dfrac{k_{y}}{k_{\parallel}}\\ -\dfrac{k_{y}k_{\parallel}}{\sqrt{\varepsilon_{1}}k_{0}k_{\parallel}}&-\dfrac {k_{x}}{k_{\parallel}}\\ \dfrac{k_{y}}{k_{\parallel}}&0\end{array} \right) \tag{4}\]  
The transmitted electric far-field of the whole beam in a \(xOy\) plane located at \(z=z_{0}>>\lambda\) is then calculated by Fourier transform of the transmitted angular spectrum. This leads to the following expression:  
\[\overrightarrow{E_{t}}(x,y,z_{0})=\dfrac{1}{4\pi^{2}}\iint\limits_{-\infty} ^{+\infty}\overrightarrow{F}_{t}(k_{x},k_{y})e^{-ik_{x}x-ik_{y}y-i(\sqrt{ \varepsilon_{2}k_{0}^{2}-k_{x}^{2}+k_{y}^{2}})z_{0}}dk_{x}dk_{y} \tag{5}\]  
Thus, injecting Equation (4) into Equation (3) and the latter into Equation (5) allows the determination of the transmitted beam, providing the knowledge of the angular spectrum of the incident beam \(\overrightarrow{F}(k_{x},k_{y})\).