A common estimation problem in signal processing assumes the following signal model \begin{equation} \mathbf{r} = \sum_{i=1}^{Q}\alpha_i\mathbf{s}\left(w_i\right)+\mathbf{n} \end{equation} where $\mathbf{s}(w)$ is a known vector function, and $\mathbf{n}\sim \mathcal{CN}(0,\sigma^2\mathbf{I})$ is white Gaussian noise with variance $\sigma^2$.

- Column vector $\mathbf{r}\in\mathbb{C}^N$ constitutes the noisy data.
- Parameters $w_1,\ldots,w_Q$ are unknown and need to be estimated.
- $\alpha_1,\ldots,\alpha_Q$ are zero-mean complex random Gaussian variables with unknown variance.
- $Q$ is the number of parameters to be estimated and is known.

A very useful family of estimators rely on some properties of the correlation matrix which is \begin{equation} \operatorname{E}\mathbf{r}\mathbf{r}^H = \mathbf{S}\left(\mathbf{w}\right) \mathbf{A} \mathbf{S}^H\left(\mathbf{w}\right) + \sigma^2\mathbf{I} \end{equation} where $\mathbf{A}=\operatorname{E}\,[\alpha_1,\ldots,\alpha_Q]^T[\alpha_1,\ldots,\alpha_Q]^*$, $\mathbf{S}(\mathbf{w})=[\mathbf{s}(w_1),\ldots,\mathbf{s}(w_Q)]$, $(\cdot)^T$ is the transpose operator, $(\cdot)^*$ is the conjugate operator, and $(\cdot)^H$ is the transpose and conjugate operator. In practice, the correlation matrix is estimated using the sample covariance matrix defined as \begin{equation} \hat{\mathbf{R}}=\frac{1}{K}\sum_{k=1}^{K}\mathbf{r}_k\mathbf{r}_k^H \end{equation} where $\mathbf{r}_k$ are independent realizations of $\mathbf{r}$.

Denote $\lambda_1,\ldots,\lambda_N$ and $\mathbf{v}_1,\ldots,\mathbf{v}_N$ the eigenvalues and eigenvectors of $\hat{\mathbf{R}}$, respectively. The eigenvalues are ordered from largest to smallest (all are non-negative because $\hat{\mathbf{R}}$ is hermitian). Then in such case it turns out that \begin{equation} \begin{bmatrix} \mathbf{v}_1 \cdots \mathbf{v}_Q \end{bmatrix} \operatorname{diag}\left(\lambda_1-\sigma^2,\ldots,\lambda_Q-\sigma^2\right) \begin{bmatrix} \mathbf{v}_1 \cdots \mathbf{v}_Q \end{bmatrix}^H \rightarrow \mathbf{S}\left(\mathbf{w}\right) \mathbf{A} \mathbf{S}^H\left(\mathbf{w}\right) \end{equation} in probability when $K\rightarrow +\infty$.

Without going into further unnecessary detail my question is

**Is it known the exact/asymptotic/approximated distribution of $\begin{bmatrix}
\mathbf{v}_1 \cdots \mathbf{v}_Q
\end{bmatrix}
\operatorname{diag}\left(\lambda_1,\ldots,\lambda_Q\right)
\begin{bmatrix}
\mathbf{v}_1 \cdots \mathbf{v}_Q
\end{bmatrix}^H$?**

As far as I know the sample covariance matrix has a Wishart distribution, but I do not know what is the distribution of such *signal* covariance matrix.