## 2 Materials and Methods
### Proposed CSK-Based UM-VLC Channel Model  
Consider a set of red, green, and blue LEDs with the same characteristics except for different SPDs, given by \(\Phi_{1}(\lambda)\), \(\Phi_{2}(\lambda)\), and \(\Phi_{3}(\lambda)\), respectively, spatially grouped together in such a way that their positions in space can be approximated to be the same from a great enough distance, and three PDs, each with spectral responsivity \(R(\lambda)\) and surface area \(A_{R}\) and filtered by either a red, green, or blue optical filter with spectral gain \(\mathcal{G}_{1}(\lambda)\), \(\mathcal{G}_{2}(\lambda)\), or \(\mathcal{G}_{3}(\lambda)\), respectively, spatially grouped so that their positions in space can be approximated to be the same from a great enough distance. If the distance \(d\) between the LEDs and the PDs is large enough, assuming that a unique transmitter position and a unique receiver position is a reasonable approximation, then the electro-optical UM-VLC DC channel gain from the \(i\)th LED to the \(j\)th PD, with \(i,j=1,2,3\), is given by  
\[H_{\text{el}}[i,j]=H_{\text{LOS}}^{\text{el}}[i,j]+H_{\text{NLOS}}^{\text{el} }[i,j]\,, \tag{91}\]  
where \(H_{\text{LOS}}^{\text{el}}[i,j]\) is the DC gain of the LOS link, which can be expressed as  
\[ \begin{split} H_{\text{LOS}}^{\text{el}}[i,j]&= \frac{1}{P_{T,i}}\int R(\lambda)\mathcal{A}_{0,j}(\lambda)\Phi_{i}(\lambda)\, d\lambda\,,\\ &=\mathcal{C}\frac{(m+1)A_{R}}{2\pi d^{2}P_{T,i}}\,\Omega_{\text {LOS}}\Psi_{\text{LOS}}\int R(\lambda)\mathcal{G}_{j}(\lambda)\xi_{0}(\lambda )\Phi_{i}(\lambda)\,d\lambda\,,\end{split}  \tag{92}\]  
where \(P_{T,i}=\int\Phi_{i}(\lambda)\,d\lambda\). Similarly, \(H_{\text{NLOS}}^{\text{el}}[i,j]\) is the DC gain of the single-hop NLOS link, which can be expressed as  
\[ \begin{split} H_{\text{NLOS}}^{\text{el}}[i,j]&= \sum_{w=1}^{W}H_{\text{NLOS, }w}^{\text{el}}[i,j]\,,\\ &=\frac{1}{P_{T,i}}\sum_{w=1}^{W}\int R(\lambda)\mathcal{A}_{w,j} (\lambda)\Phi_{i}(\lambda)\,d\lambda\,,\\ &=\frac{(m+1)A_{R}}{2\pi P_{T,i}}\\ &\quad\times\sum_{w=1}^{W}\mathcal{C}_{w}\frac{A_{\text{real},w} }{d_{1,w}^{2}d_{2,w}^{2}}\,\Omega_{\text{NLOS},\,i}^{(w)}\Psi_{\text{NLOS}}^{ \text{(w)}}\int R(\lambda)\,\rho_{w}(\lambda)\mathcal{G}_{w,j}(\lambda)\xi_{w }(\lambda)\Phi_{i}(\lambda)\,d\lambda\,.\end{split}  \tag{93}\]  
Note that if we consider that optical filters do not depend on incident angle \(\theta_{w}\), then we can relax the definition to \(\mathcal{G}_{w,j}(\lambda)\equiv\mathcal{G}_{j}(\lambda)\). Furthermore, \(\Omega_{\text{NLOS},\,i}^{(w)}\) depends on the \(i\)th LED, since we consider that the light that bounces off the \(w\)th reflector will have different angles of irradiance depending on the LED from which the light comes. This is to account for the irregularity of the underground tunnel walls, where a slight displacement of the light source can have great effects on the angle of reflection off the wall. Since the LEDsare slightly displaced from each other, it is safe to assume that light from different LEDs will bounce off the walls independently. In practice, we model the \(w\)th reflector differently for each LED, such that it has a rotation angle \(\alpha_{1,w}\) and a tilt angle \(\beta_{1,w}\) for light coming from the red LED; a rotation angle \(\alpha_{2,w}\) and a tilt angle \(\beta_{2,w}\) for light coming from the green LED; and a rotation angle \(\alpha_{3,w}\) and a tilt angle \(\beta_{3,w}\) for light coming from the blue LED. For each \(w\)th reflector, these six angles are modeled as uniform random variables, drawn independently.