Let $\mathbf{S}^n$ denote the space of symmetric, real-valued $n \times n$ matrices.
Consider the closed convex set
$$ \mathcal{C} := \{(X, x) \in \mathbf{S}^n \times \mathbf{R}^n : X \succeq xx^T, ~ \mathbf{tr}(X) \leq 1\},$$ where $\succeq$ above denotes the positive semidefinite (Lowner) order, and $\mathbf{tr}(\cdot)$ denotes trace.
I would like to compute the Euclidean projection onto this set, i.e., I wonder if there is a closed form for the following operator, $\mathrm{proj}: \mathbf{S}^n \times \mathbf{R}^n \to \mathcal{C}$, which is variationally given by $$ \mathrm{proj}(Z, z) = {\mathrm{argmin}}_{(X, x) \in \mathcal{C}} \left(\frac{1}{2}\|Z - X\|_F^2 + \frac{1}{2}\|z - x\|_2^2 \right), $$ where above $\|\cdot\|_F$ denotes the Frobenius norm.
First note since $ X \succeq xx^T \Leftrightarrow \begin{bmatrix} X & x \\ x^T & 1 \end{bmatrix} \succeq 0 $, it makes sense to look at the problem as one in $S^{n+1}$. Let $P_1, P_2$ be linear projections from $S^{n+1}$ to $S^{n} \times R^n$ and $R$ respectively: $$ P_1 \left( \begin{bmatrix} X & x \\ x^T & y \end{bmatrix} \right) := (X,x), \quad P_2 \left( \begin{bmatrix} X & x \\ x^T & y \end{bmatrix} \right) := y $$
So then define the subsets of $S^{n+1}$: $$ \begin{aligned} D_1 &:= \{ W \in S^{n+1} : W \succeq 0, ~ tr(W) \le 2\}\\ D_2 &:= \{ W \in S^{n+1} : P_2(W) = 1\} \end{aligned} $$ Then the intersection of these is equivalent to $C$ in the following sense: $$ (X,x) \in C \Leftrightarrow P_1^T(X,x) + P_2^T(1) \in D_1 \cap D_2 $$ Note that it is possible to project onto $D_1$ using eigenvalue decomposition and projection onto $D_2$ is trivial, since it's just a plane. So we can almost get a closed form expression for the projection onto the intersection.
Using the notation $ \Pi_A (z) = \mathrm{argmin}_{x \in A} \|x-z\|^2 $ for nonlinear projection onto a set, we have:
$$ \Pi_C (Z,z) = P_1\left( \Pi_{D_1} \left( \begin{bmatrix} Z & z \\ z^T & y^* \end{bmatrix} \right) \right) $$ Where $y^*$ is the value satisfies the following: $$ 1 = P_2\left( \Pi_{D_1} \left( \begin{bmatrix} Z & z \\ z^T & y^* \end{bmatrix} \right) \right) $$ This can be solved numerically via a 1-dimensional root-finding method. I'm not sure if it's possible to avoid calculating a full eigenvalue decomposition every time you project onto to $D_1$.