}_{j=1}^N {{\Bar{\Bar{\mathbf{G}}}}}({\mathbf{r}}, {\mathbf{r}}_j, \omega)\cdot {\mathbf{d}}_j(\omega).

\end{array}$$ ]{} Here, the dyadic Green’s tensor ${{\Bar{\Bar{\mathbf{G}}}}}({\mathbf{r}}, {\mathbf{r}}_j, \omega)$ encodes the field at position ${\mathbf{r}}$, produced by an oscillating dipole at ${\mathbf{r}}_j$, and in vacuum is given by [@Novotny2009PrinciplesNano-optics] [ $$\begin{array}{c}

\label{eq:FTMS_Green_tensor}
{{\Bar{\Bar{\mathbf{G}}}}}({\mathbf{r}},{\mathbf{r}}',\omega) =
k \dfrac{e^{i \rho}}{4\pi }\left[\left(\dfrac{1}{\rho}+\dfrac{i}{\rho^2}-\dfrac{1}{\rho^3}\right)\mathbb I
+\left(-\dfrac{1}{\rho}-\dfrac{3i}{\rho^2}+\dfrac{3}{\rho^3}\right)\dfrac{\boldsymbol{\rho} \otimes \boldsymbol{\rho}  }{\rho^2}\right],

\end{array}$$ ]{} with $\rho \equiv |\boldsymbol\rho|\equiv k| ({\mathbf{r}} - {\mathbf{r}}')|$ and $k=\omega/c$. Note that ${{\Bar{\Bar{\mathbf{G}}}}}({\mathbf{r}},{\mathbf{r}}',\omega)$ contains both non-radiative, near-field ($\sim 1/\rho^3$) and radiative, far-field ($\sim 1/\rho$) terms.

Then, the induced dipole moment of atom $i$ is given by [ $$\begin{array}{c}

\label{eq:FTMS_coupled_dipoles}
 d_i(\omega)=\alpha_0(\omega) \epsilon_0\left[ E_{\text{in}}({\mathbf{r}}_i,\omega)+\mu_0 \omega^2 {
\displaystyle \sum 
}_{j\neq i}^{N-1} {\hat{\mathbf{x}}} \cdot {{\Bar{\Bar{\mathbf{G}}}}}({\mathbf{r}}_i,{\mathbf{r}}_j, \omega)\cdot {\hat{\mathbf{x}}}\;  d_j(\omega)\right],
\end{array}$$ ]{} where the parameter $\alpha_0(\omega)$ defines the polarizability of a single dipole. Although [Eq. \[eq:FTMS\_output\_field\]]{} and [Eq. \[eq:FTMS\_coupled\_dipoles\]]{} can describe any system of linearly-polarizable point-like dipoles [@Novotny2009PrinciplesNano-optics], e.g. dielectric nano-particles [@GarciaDeAbajo2007CollectiveMatter], in our case we focus on the response of non-absorbing, purely radiative atoms, whose resonant cross section $\sigma_{\text{sc}}=3\lambda_0^2/(2\pi)$ is the maximum set by the unitarity limit [@DeVries1998PointWaves]. In this context, the atomic polarizability reads $\alpha_0(\omega)=-3\pi /[(\Delta + i/2)k_0^3]$, where $k_0=2\pi/\lambda_0$ denotes the resonant wavevector, while $\Delta\equiv(\omega - \omega_0)/\Gamma_0$ represents the dimensionless detuning between the input beam frequency $\omega$ and the atomic resonance $\omega_0$.

