Rotation sensing at the ultimate limit

Conventional classical sensors are approaching their maximum sensitivity levels in many areas. Yet these levels still are far from the ultimate limits dictated by quantum mechanics. Quantum sensors promise a substantial step ahead by taking advantage of the salient sensitivity of quantum states to the environment. Here, we focus on sensing rotations, a topic of broad application. By resorting to the basic tools of estimation theory, we derive states that achieve the ultimate sensitivities in estimating both the orientation of an unknown rotation axis and the angle rotated about it. The critical enhancement obtained with these optimal states should make of them an indispensable ingredient in the next generation of rotation sensors that is now blossoming.


Introduction
Precise rotation sensing is a critical requisite for applications as diverse as inertial navigation [1,2], geophysical studies [3], and tests of general relativity [4,5,6], to cite but a few.
Since the first gyroscope, invented in 1852 by Foucault as an extension of his famous pendulum used to demonstrate the Earth's rotation, rotation sensing has experienced a formidable transformation. Depending on the application and its environmental conditions, various types of sensors have been reported. Some of the most commonly used are based on Sagnac interferometers [7], although they suffer from limited sensitivity and stability. Matter-wave interferometers have demonstrated superior performance [8,9,10], but at the price of a substantial increase in complexity. Trapped atoms in a guiding potential have also been proposed to improve rotation precision [11,12].
Irrespective of the technique, sensors are reaching performance levels where quantum effects come into play [13] and, accordingly, should be re-examined from a full quantum perspective: this is the main goal of this paper. Actually, quantum metrology ascertains the ultimate bounds on the achievable measurement precision and identifies states that would be optimal for those measurements [14]. The main tool for those tasks is the quantum Cramér-Rao bound, which bounds the covariance matrix between parameters being estimated by the inverse of the quantum Fisher information matrix [15,16].
When the axis of the rotation is known, the only endeavour is to estimate the angle of rotation. This is entirely equivalent to determining a phase in interferometry [17], a paradigmatic example of single-parameter estimation. This topic is well understood [18,19]: classical states have a precision limited by shot noise, leading to uncertainties scaling as 1/ √ N , where N is the number of particles involved in the measurement [20,21]. However, one can find special quantum states saturating the quantum Cramér-Rao bound: they achieve the so-called Heisenberg limit, in which measurement uncertainties scale as 1/N [22,23,24,25,26,27].
In general, however, a rotation is characterized by three parameters [28]: either the two angular coordinates of the rotation axis and the angle rotated around that axis, or the Euler angles [29]. We thus face the problem of simultaneously estimating of multiple parameters, which has been considered as a distinguished feature combining classical and quantum aspects of uncertainty: the noncommutativity of quantum theory leads to nontrivial tradeoffs that are absent in classical and in singleparameter estimation problems. This observation led to a development of new multiparameter bounds that constitute the backbone of a new and lusty research line [30,31,32,33,34,35,36].
In this paper, we discuss optimal states for the simultaneous estimation of the three parameters of a rotation. These states attain the Heisenberg limit, confirming that they always achieve a significant enhancement in measurement precision relative not just to the classical case but also to the best single parameter schemes. Finally, we briefly discuss their practical implementation.

Classical estimation
Estimation theory deals with devising schemes that extract as precisely as possible the value of an unknown parameter and therefore lies in the realm of metrology. From a physical perspective, a typical estimation process can be roughly divided into three stages, which are schematized in figure 1: probe preparation, interaction with the system, and probe readout.
The measured data are always affected by noise, so they are effectively represented by a stochastic variable. In some cases, the experiment cannot be modeled mathematically and the use of nonparametric estimation is necessary [37]. However, it is always more efficient to find a convenient model for the studied experiment. Such a parametric estimation will be the main focus of this paper.
In classical (or frequentist) estimation, the unknown parameter is taken to be deterministic and constant during the experiment: randomness is thus solely due to noise. In contradistinction, Bayesian estimation assumes that the unknown parameter is itself a random variable distributed according to some prior probability distribution. Whenever this distribution is unknown, it is safer to perform classical estimation, which will be our primary goal.

Initial state
Final state Dynamical process Measurement Estimator Figure 1. Scheme of a parameter estimation protocol. A probe is prepared, then it interacts with a dynamical process that imprints the unknown parameter θ onto the probe. The resulting state, encoding the information about θ, is measured and generates outcomes x. Based on the outcomes x, a suitable estimator θ provides an estimate of the parameter θ.

Single-parameter estimation
Let us first consider the problem of estimating the value of a single parameter θ from the observed measurement result x, where the latter is taken to be a continuous random variable (the discrete case can be dealt with much in the same way). In general, experiments consist of N trials yielding a set of outcomes x = (x 1 , . . . , x N ) , where the superscript denotes the transpose. The inference of the parameter is related to the measurement outcomes through some conditional probability density that is dictated by a model of the process and that we denote by P (x|θ). The N random variables are typically assumed to be independent and identically distributed (iid). Identically distributed implies that each of these N random variables is governed by the same parameter θ, whereas independence implies that the joint distribution of all these N random variables is given by a product: (2.1) To extract the parameter θ from the data x we use an estimator θ(x), which is a function of the observed data only. This estimator is said to be unbiased if its mean value coincides with the true value of the unknown parameter (i.e., there is no systematic error in the estimation) where we denote by the expectation value of the estimator θ(x) with respect to the probability distribution P (x|θ). A vast number of estimators have been formulated by proposing ad hoc functions θ(x). The reader is referred to the excellent textbooks on parametric estimation [38,39,40]. One of the most widely adopted is the maximum likelihood (ML) estimator, defined as Here, the likelihood P (θ|x) must be interpreted as the probability of θ being the true value given the data set x. This estimator is unbiased in the asymptotic limit of many independent samples. Actually, it just picks the value of θ that makes x the most likely observation; in other words, it bets on the winner. A natural figure of merit to quantify the performance of an estimator is its mean square error (MSE) (2.5) This definition depends on the unknown parameter θ, but the MSE is a priori a property of the estimator. Actually, it can be decomposed as 6) where the variance and the bias are, respectively, Bias 2 The minimum MSE estimator finds a trade-off between the bias and the variance for every value of θ. Unfortunately, the bias is often a function of θ and, consequently, the minimum MSE estimator is generally not realizable. In general, any estimator depending on the bias will be unfeasible. This limitation prompts one to focus uniquely on unbiased estimators, with the resulting estimator usually referred to as the minimum variance unbiased (MVU) estimator. There are procedures for finding the MVU estimator [40]. Unfortunately, they are often tedious and sometimes fail to produce the MVU estimator, whose existence is not even guaranteed. Fortunately, the ML estimator is known to approximately provide the MVU estimator under mild regularity conditions.
A fundamental tool to characterize the achievable bounds on estimation uncertainties is the Fisher Information (FI) [41], which captures the amount of information encoded in the output probabilities. For the random variable x, it is defined by where the alternative form in the second line can be obtained by a direct computation. We keep the subscript x as a reminder that the FI hinges upon the data x, which in fact depends on both the probe state and the measurement. The choice of the function ln P (x|θ) is meaningful because it is additive for independent samples due to the logarithm. Intuitively, the FI quantifies the sensitivity of a system to a change in θ: a larger amount of information is associated with larger variations of the output probabilities. This intuition was formalized with a time-honoured result called the Cramér-Rao bound (CRB) [42,43]. It links F x (θ) with the ultimate bound achievable by the variance of any arbitrary unbiased estimator, with N identical and independent probes and measurements: . (2.9) The CRB applies only to well behaved probability distributions that satisfy the following regularity condition [40] An estimator that saturates the inequality (2.9) is said to be efficient, which can be guaranteed by the condition (2.11) In the limit of a large number of measurements, the ML estimator is efficient, for its distribution normally converges to the real value with a variance that saturates the CRB.
The derivation of the CRB relies on the Cauchy-Schwarz inequality using ∂ ln P (x|θ)/∂θ, which is called the score function, as the basic object. Other choices of the score give different bounds, such as the Hammersley-Chapman-Robbins [44,45], Bhattacharya [46], and Barankin [47] bounds.
To conclude, we point out that in assessing the performances of an estimator, a relevant factor is the scaling of the variance with the mean value rather than the absolute value of the variance. This feature may be quantified by means of the signalto-noise ratio [30] (for a single detection event) which is larger for better estimators. Using the CRB, one easily derives that the signal-to-noise ratio of any estimator is bounded by the quantity R The parameter θ is effectively estimable when the corresponding R x ( θ) is large.

Multiparameter estimation
The central task is now estimating a series of D unknown parameters θ = (θ 1 , . . . , θ D ) , which are assessed through a set of estimators θ = ( θ 1 , . . . , θ D ) , after the measurement x = (x 1 , . . . , x D ) . As before, we assume N trials that give N distinct iid random variables x. Each parameter θ i , with i = 1, . . . , D, represent a physical quantity. The figure of merit is now the covariance matrix, with elements (i, j = 1, . . . , D) which quantifies the sensitivity relative to each parameter, while also taking into account the possible correlations between the parameters. Much in the same way, the FI is generalized to the FI matrix F x (θ) with elements which is symmetric and positive. This matrix can be interpreted as a metric on the parameter space and transforms as [33] for Jacobian matrix J ij (θ, ϑ) = ∂θ i /∂ϑ j governing the transformations between parametrizations.
For an unbiased estimator, the CRB in the multiparameter case is generalized to the following matrix inequality An important example is the Gaussian (normal) distribution, which allows for a closed form of F. If the probability distribution is characterized by mean values µ and covariance matrix Σ, the FI matrix reads

Quantum estimation
We next review the quantum setting of the problem. The general scheme is illustrated at the bottom of figure 1. From a quantum perspective, this protocol involves a twostep optimization problem: one must first make a smart choice of the probe state that is sensitive to the parameter, and then pick an appropriate measurement that maximizes the information extracted from the probe.

Single-parameter estimation
We represent the probe state by a density operator . The single parameter θ is encoded via a quantum channel E θ , whose action on the state is the transformation [48]. For simplicity, we will restrict our attention to unitary channels, so that θ = U θ U † θ , with the unitary operator U θ expressed as where G is a selfadjoint operator that is called the generator of the transformation. In order to estimate θ, we perform a measurement that is represented by some positive operator-valued measure (POVM) {Π x }. The latter comprise a set of positive semidefinite, selfadjoint operators that resolve the identity [49]; that is, By performing this measurement, we obtain a statistical distribution that, according to Born's rule, is given by Afterward, what remains is to obtain the best estimate of θ given P (x|θ). This can be accomplished with the basic tools outlined in the previous section, albeit one has to guarantee the positivity of the quantum state. In other words, quantum estimation can be seen just as classical estimation supplemented with the constraints imposed by positivity.
There is an infinite number of possible POVMs that we can consider. It is therefore natural to ask whether there is an optimal measurement that should be performed on θ . The quantum Fisher information (QFI) is defined precisely for this purpose [50]: which, as indicated by its subscript, depends exclusively on the initial probe state . By its very same definition, we have Q (θ) ≥ F x (θ) and so we have a quantum Cramér-Rao bound (QCRB) given by .
Here, the variance must be computed using the probability density (3.3) associated with Born's rule. The right hand side of (3.5) thus represents the ultimate achievable precision regardless of the measurement. Helstrom, using fairly elementary arguments [19], showed an explicit way to compute Q (θ); it reads where the variance of an operator is Var (L) = Tr( L 2 )−[Tr( L)] 2 , L θ is the so-called symmetric logarithmic derivative, defined implicitly via and {·, ·} stands for the anticommutator {A, B} = AB + BA. The QFI has a number of interesting properties. First, it is convex in the quantum states: given any two states 1 and 2 , we have with p 1 + p 2 = 1. Moreover, the QFI is additive for independent measurements; that is, Second, for unitary evolutions, the QFI does not depend on the position along the orbit of U θ ; i.e., it is the same for the state as for θ = U θ U † θ . In such a case, the QFI simplifies to [30] where p i and |i are the eigenvalues and eigenvectors of , respectively. For pure initial states, = |ψ ψ|, under unitary evolution, a simpler expression for the QFI is where and g max and g min are the maximum and minimum eigenvalues of G, respectively. The QFI is thus maximal for pure states that maximise the variance of G.
The QFI is intimately connected with the distinguishability of a probe for small variations of the parameter [51]. The distinguishability between two states, 1 and 2 , can be quantified by the normalized Bures distance [52] [53]. Then, given two states θ and θ+δθ , obtained by an infinitesimal change in the parameter, one has [33] except for pointwise differences when θ changes rank [54,55]. Therefore, the more distinguishable a state is, the greater the QFI and the sensitivity of the state to the parameter θ.
A final point of paramount importance is to determine measurements that reach the ultimate precision and thus saturate the QCRB. This is equivalent to finding a POVM such that the associated FI is equal to the corresponding QFI for the probe state. In the single parameter case, it is always possible to saturate the QCRB by simply taking projectors onto the eigenstates of the symmetric logarithmic derivative L θ [19].

Multiparameter estimation
The extension to multiparameter quantum estimation looks superficially similar to the classical case. The QFI is generalized to the real-valued symmetric QFI matrix with components where L i is the symmetric logarithmic derivative with respect to the parameter θ i . For the particular case of pure states, this reduces to where |∂ θi ψ θ = ∂|ψ θ /∂θ i . This QFI matrix verifies convexity and additivity properties [56]. It follows that a multiparameter QCRB can then be formulated as Summing over the diagonal elements of this matrix inequality we arrive at For the same resources, it is now established that the simultaneous quantum estimation of multiple parameters provides better precision than estimating them individually [57]. Unfortunately, the possibility of attaining the ultimate quantum bounds for the simultaneous estimation is not guaranteed [58,59,60,61]: the corresponding optimal measurements may not commute, thus making their implementation impossible.
The topic of multiparameter quantum estimation is nowadays quite an active field of research that has emerged as the confluence of several disparate yet interconected developments in different fields. The interested reader should consult recent reviews on the topic [30,31,32,33,34,35,36]. For our purposes here, it is enough to mention that a sufficient condition for the joint estimation is that the operators L i commute. A weaker condition is provided by the following constraint [60] For pure states, there exists a necessary and sufficient condition for the saturation of the QCRB: if Q (θ) is invertible, the QCRB can be saturated if and only if When the condition (3.17) is not met, there exists a tighter bound based on the socalled right logarithmic derivative [62]. However, this operator is not directly linked to any measurement, in contradistinction to the L i .

Rotations and the Majorana stellar representation
Let us consider a rotation defined by its axis n = (sin Θ cos Φ, sin Θ sin Φ, cos Θ) and angle θ ∈ [0, π), measured according to the right-hand rule. We will use the compact notation θ = (θ, n) to denote these parameters. Under this rotation, a point of the system with coordinates r rotates to a new position r , given by Henceforth, the latin indices {i, j} will run over the coordinate indices {1, 2, 3}, ε ijk is the totally antisymmetric unit vector, and summation over repeated indices is assumed.
If one particle is in the state |r , representing its localization at the point r, after a rotation R θ its state must be proportional to |R θ r . A celebrated theorem of Wigner ensures that this action is implementable by unitary or antiunitary operators [63]. The latter possibility can be discarded because each rotation is the square of some other rotation. In consequence, where R θ is the unitary representation of the rotation (defined up to a global phase). Therefore, we have that which shows that the wave function for a particle without internal degrees of freedom transforms under rotations as a scalar. Equation (4.4) allows us to find the explicit unitary operator R θ . A direct calculation shows that [64] which is defined up to a ± sign due to the projective character of the mapping R → R that is globally unavoidable [65]. Here, the operator J is the generator of all the effects of a rotation and is the total angular momentum of the system. It is an observable composed of the total orbital angular momentum L and the total spin S. The three components {J i } satisfy the commutation relations of the Lie algebra su(2) The irreducible representations are labelled by J and spanned by the states {|Jm }, which are the simultaneous eigenstates of J 2 and J 3 (with = 1): Since m = −J, . . . , +J, these states span a (2J + 1)-dimensional Hilbert space we denote by H J . In what follows, we will assume that we work in H J , unless explicitly stated otherwise. The notion of a Majorana representation [66] will prove to be extremely convenient for our purposes. In this representation, a pure state corresponds to a configuration of points on a sphere, a picture that makes a high dimensional Hilbert space easier to comprehend. The idea can be presented in a variety of ways [67,68], but the most direct one is, perhaps, by realizing that every pure state |ψ ∈ H J can be mapped onto the polynomial where ψ m = Jm|ψ are the amplitudes of the state in the angular momentum basis. Up to a global factor, |ψ is determined by the set {z i } of the 2J complex zeros of ψ(z) (suitably completed by points at infinity if the degree of ψ(z) is less than 2J). A nice geometrical representation of |ψ by 2J points on the unit sphere (often called the constellation) is obtained by an inverse stereographic map Notice that the location of the stars has an operational meaning: a spin system, say, cannot be measured to have spin up along a direction that points away from a star on the sphere. To illustrate how this representation works in practice, we will examine a few relevant examples. The first one is that of SU(2) or Bloch coherent states. They can be defined, among other equivalent ways, as [69,70] where J ± = J 1 ± iJ 2 are ladder operators, and n denotes a unit vector in the direction of spherical angles (θ, φ) defined in (4.9). Coherent states are precisely eigenstates of the operator J · n. For the state |z 0 we have so the polynomial has a single zero at z = −z 0 with multiplicity 2J. Consequently, the constellation collapses in this case to a single point diametrically opposed to n 0 . The Majorana constellations for the angular momentum basis |Jm can be easily inferred from the polynomial so they consist of J ± m stars at the north and south poles, respectively. Another relevant set of states are the so-called NOON states, defined as [71] with associated polynomials (4.14) The zeros are the 2J roots of unity, so the Majorana constellations have 2J stars placed around the equator of the unit sphere with equal angular separation between each star. We can generalize the NOON states to Bloch cat states [72] They interpolate from NOON states (|z| = 1) to angular momentum eigenstates |J ± J or |J ± (J − 1) (|z| = 0, ∞). The Majorana polynomial is Since the most classical states have the most concentrated constellation, one might naively think that the most quantum states have their 2J stars distributed in the most symmetric fashion on the unit sphere, and this is the caase. This constitutes the realm of the so-called Kings of Quantumness [73,74], initially dubbed anticoherent states [75]. In a sense they are the opposite of Bloch coherent states: while the latter correspond as nearly as possible to a classical spin vector pointing in a given direction, the former point nowhere; i.e., the average angular momentum vanishes and the fluctuations up to given order M are isotropic [76,77,78]. Their symmetrical Majorana constellations herald their isotropic angular momentum properties; these states are the most sensitive for rotation measurements. A few examples of these constellations are shown in Figure 1.

Known rotation axis
We first consider estimating the angle θ by which a system is rotated about a known axis n. Since the angular momentum operators generate rotations via we can always take n as being directed along the z axis. The generator is then G = J 3 , with maximal and minimal eigenvalues corresponding to eigenstates |JJ and |J − J .
The optimal probe states are those maximizing the variance of J 3 . This requires an equal superposition of the eigenstates |JJ and |J − J ; that is, the NOON states treated in the previous section. Estimating a rotation angle for a known rotation axis is mathematically equivalent to phase estimation, explaining the ubiquity of NOON states in that context, with Heisenberg-scaling performances on average.
When the rotation axis is an arbitrary n, the optimal probe states are the rotated versions of NOON states; i.e., (|n +|−n )/ √ 2, as the rotated version of the north pole |JJ is precisely a Bloch coherent state. Instead of producing a different state for each axis, one may be interested in a probe state that has the best average performance for any n. There are a few ways of quantifying such performance. First, one may consider states whose average Fisher information is largest: Var θ (J i ).

(5.2)
Optimizing this quantity is the same as optimizing the signal-to-noise ratio for a given total angular momentum J. The minimum is achieved by states whose angular momentum projection is maximal in any fixed direction; i.e., they satisfy Tr( θ J) = Jn . For pure states, they are |n andQ min = J. The maximum is achieved by states whose angular momentum projection vanishes in all directions; i.e., Tr( θ J) = 0, withQ max = 4 3 J(J + 1). In this sense, the most sensitive probe states are those whose classical angular momentum features are hidden (they include the NOON states).
One can also find states that optimize the average variance of the estimated parameter by averaging over the QCRB (3.5): .
In that case, the average over all angles diverges for the states R θ |Jm for any θ , ‡ while it is minimized by states satisfying with value 3 4 J(J + 1). The states achieving the minimal average variance are those whose first and second order angular momentum properties are isotropic, which ‡ Incidentally, this supports the conjecture in Ref. [79] that the average over all angles diverges for more than just the states |n , but we have not proven the nonexistence of other such states.
are known as Kings of Quantumness. These have Heisenberg-scaling uncertainties on average. NOON states, as an intermediate, achieve an average variance of arctan( , which scale more poorly with J than the Kings but better than the states |n . This shows how the problem of rotation estimation requires more nuance than that of phase estimation [79]. The outstanding performance of estimating rotation angles about various axes has been demonstrated using light's orbital angular momentum degree of freedom in Ref. [80].
It is worth mentioning that the Kings of Quantumness also have isotropic higherorder moments. Since the only criterion generally accepted for assessing the optimality of an estimator is whether its variance saturates the QCRB, the standard approach might overlook relevant information that may be present in the complete parameter distribution, and that can be captured by looking at higher-order moments [81].

Unknown rotation axis
The rotation axis n is not always known a priori. Then, we must consider the multiparameter estimation problem of determining all three components of θ. We mainly focus on our θ parametrization due to its physical importance, but our results are applicable to the various parametrizations one can find in the literature. Of note, one can always find the QFI matrix for a different parametrization by using equation (2.15).
Before determining and optimizing the QFI matrix for this multiparameter scenario, we mention that an alternative approach to optimizing states for probing rotations about arbitrary axes is that of optimal quantum rotosensors [82,83]. These are probe states that are optimally sensitive to determining whether a system has undergone a rotation by known angle θ about unknown axis n [83]. This is done by finding states that are, on average, the most changed by a rotation for a fixed angle. Then, the most sensitive states for identifying the presence of a small rotation angle are the Kings, of a large rotation angle are the states |n , and of intermediate rotation angles are other states. Optimal quantum rotosensors can be used to determine whether or not a rotation is present; if present, one can use the results presented here to optimally determine the parameters of the rotation.
The three rotation parameters being estimated correspond to three generators defined as G k = (i∂ θ k R θ )R † θ . Computing these generators is nuanced because [∂ θ k (J · n), J · n] = 0. To fix the problem, we express the generators as [84] where ω = θn, and the final expression is guaranteed to be a linear combination of angular momentum operators by the Baker-Campbell-Hausdorff formula [85] or equivalently by the properties of rotation operators. This allows us to compute the vectors g k (θ) using any representation of angular momentum, for example by computing the vectors using Pauli matrices as done in Ref. [86]. We find the vectors g k (θ) to be given by Then, the QFI matrix can be expressed as G(θ) = (g θ g Θ g Φ ) encompasses the angular information of the three parameters and the sensitivity characteristics of the rotated state. Since the covariances may vary with θ, we can consider the transformation property and incorporate the rotation matrices into the angular matrices withG(θ) = R θ G(θ). Equation (5.10) has the same form as the change-ofparametrization formula for the QFI matrix in equation (2.15). One can easily switch to a different parametrization: the Jacobian governing this change will, similar tõ G(θ), only depend on geometric properties of the coordinate systems, and so the final QFI matrix will always have the form Q ψ (ϑ) = 4G (ϑ) C(ψ)G(ϑ) for any set of parameters ϑ. The crucial point is that in order to find states optimally suited for estimating arbitrary unknown rotations one must optimize the sensitivity covariance matrix C(ψ). First, we identify when the QFI matrix can and cannot be inverted. A singular QFI matrix implies that the triad of parameters cannot be estimated for a given state and a given parametrization (see Appendix A for further discussion of singular QFI marices). This is the case when either C(ψ) is singular, implying that the state will never be useful for estimating all of the parameters of a rotation, or when G is singular, implying that the coordinate system is singular at that specific set of parameters regardless of the probe state. As discussed in Ref. [87], singularities in one coordinate system can be alleviated for specific parameters by switching to a new coordinate system. For example, the spherical and ZY Z Euler angle parametrizations are singular for a rotation by 0, but the Cartesian and XY Z Euler angle parametrizations are invertible there, and the Cartesian parametrization is singular for rotations by 2π and the XY Z Euler angle parametrization when the Y rotation is by ±π/2.
It is straightforward to show that C(ψ) is singular if and only if the probe state is an eigenstate of some angular momentum projection; that is, proportional to R θ |Jm . Similar to the single-parameter scenario, where eigenstates of the angular momentum projection operators were the least useful for estimating rotations about unknown axes, we see that states with any definite angular momentum projection cannot be used for simultaneously estimating all three parameters of a rotation.
To find the most sensitive states we seek to minimize the trace of the inverse of C(ψ). This is straightforward to optimize because, for any symmetric, invertible with the trace of the inverse achieving the minimum only for the Kings of Quantumness. We again see that having isotropic angular momentum up until second order makes a state most sensitive to arbitrary rotations about arbitrary axes. Next we showcase how to use these states to saturate the QCRB for a particular parametrization.

Optimal measurements
Saturability of the QCRB hinges upon the expectation values of the commutators of the generators defined by equation (5.5). Fortunately, for states with isotropic covariance matrices, the expectation values of the commutators are guaranteed to vanish. In fact, these expectation values will vanish for all states whose average angular momentum vanishes, in any parametrization with generators defined by G k = J · g k . The Kings thus guarantee that all three parameters can be simultaneously estimated at a precision saturating the QCRB for any triad of rotation parameters. We can find a projection-valued measure (PVM) that saturates the QCRB by considering projections onto the three states proportional to G k |ψ θ and onto the original rotated state |ψ θ . Because the vectors g k are orthogonal, these four states are all mutually orthogonal for the Kings. A PVM consisting of the state |ψ θ and any orthonormal combination of the states G k |ψ θ suffices to produce a classical FI matrix that equals the QFI matrix. When the probe states only have their first-order angular momentum projections vanish, one must use a Gram-Schmidt orthogonalization procedure among the states G k |ψ θ [61].
The QCRB can be saturated in the asymptotic limit, when θ tends to the true value θ. In this case, one has a good estimate |ψ θ that differs from the true state by some small, unknown rotation angle ε 1 about some unknown axis u: Defining four orthonormal projectors Π k = |φ k φ k | by the states where we can choose the three vectors v k to be orthonormal, we find the expectation values of the projectors to be (5.14) In this limit, that the probability p 0 approaches unity and the probabilities p k are vanishing shows that the estimates θ are approaching the true values θ. In addition, we have: which holds for any state whose angular momentum projection vanishes and the final equality is valid for any state whose second-order angular momentum features are isotropic. From this we can immediately compute (the derivatives of p 0 vanishing to lowest order in ε) where G ij (θ) = v k · g j . A PVM comprised of, say, the state |ψ θ and J k |ψ θ for k ∈ (1, 2, 3) will thus saturate the QCRB for optimal states.
Thus far, we have been considering only ideal scenarios. The fundamental bounds change in the presence of losses and decoherence. Actually, NOON states are extremely susceptible to losses: in the photonic case, a loss of, e.g., a single photon makes the state completely useless for sensing. The Kings of Quantumness behave much better in the presence of losses.
In practice, however, those optimal states are notoriously hard to prepare, apart from the regime of very small J. For large photon numbers, the only practically accessible states of light are squeezed states [88]. We will see in Sec. 6 why these states are indeed optimal in the asymptotic limit of very large J.

Suboptimal measurements
The PVM described in the above section may be challenging to implement experimentally. Easier is to project the rotated state onto a set of Bloch coherent states for various directions and to reconstruct the rotation parameters from these measurements.
The set of projections constitute the Husimi function [89]. Knowledge of all of the projections q n is tantamount to knowledge of the rotated state |ψ θ , but such information is redundant: it suffices to sample the function at a few locations q n and use these results to orient the Husimi function and thus estimate the rotation parameters. This method can be used for any probe state and amounts to a positive operator-valued measure (POVM) rather than a PVM because the states |n are not mutually orthogonal for differing n.
At how many locations must the Husimi function be sampled to uniquely orient it? In general, the answer depends on the probe state in question and the locations being sampled, but we argue using concepts applied from geographical positioning systems that four projections q n should suffice for this orientation problem.
Projecting the rotated state |ψ θ onto an arbitrary Bloch coherent state |n 1 amounts to sampling the Husimi function at n 1 . The value of this first projection q n1 defines a set of level curves, and the state must be oriented in such a way that n 1 lies on one of these curves. Rotating the state along any of these level curves will produce the same value q n1 .
Projecting the rotated state onto another coherent state |n 2 defines another set of level curves. Rotating the state along these level curves again produces the same value q n2 , so in general we expect there to be multiple intersection points for orienting the Husimi function such that n 1 lies along a curve q n1 and n 2 lies along a curve q n1 .
Barring pathological cases, a third projection uniquely specifies one of the above intersection points for orienting the Husimi function. A fourth projection helps deal with pathological cases, and a fifth projection helps with normalization. Viewed in a different way, the set of d angular coordinates n i can be rigidly rotated until the d projections q ni match the given state |ψ θ . The obvious pathological cases are those for which the Husimi functions at a set of d angular coordinates are equal to the Husimi functions at a rigid rotation of those angular coordinates; this is easily mitigated by knowing the nonoriented Husimi function ahead of time and thus knowing its pathologies.
Any set of projections that avoids pathological cases will suffice for discerning the rotation parameters. Uncertainty in the values {q n } leads to "fuzziness" that broadens the level curves for each projection, so increasing the number of projections increases the overall precision, just like increasing the number of satellites measuring distances increases precision for global positioning systems. The limit of measuring projections in all directions is tantamount to performing full quantum state tomography on the rotated probe state |ψ θ . As an example, consider the state in Figure 3, which has the randomly-chosen components {0.000394688 + 0.409134i, 0.0324599 + 0.0448131i, 0.494021 + 0.484609i, 0.483644 + 0.114779i, 0.100279 + 0.305783i} in the |2m basis. The figure depicts a scaled version of Husimi function for all values of n. We can consider various level curves of this function, with their uncertainties, by plotting Gaussians centered at particular values q n . For a certain rotation R θ |ψ , three separate measurements q ni at three values of n i , the white semicircles in Figure 4, produce values q n1 = 1/10, q n2 = 3/10, and q n3 = 5/10. This tells us that the Husimi function must be oriented in such a way that the three directions lie on the three level curves.
The projections {q n } can be used in, for example, maximum likelihood estimation schemes to determine the angular parameters θ. Actually, we can sketch the usefulness of this POVM in the language of Fisher information. To that end, we assume that the sampled directions are enough to provide a (approximate) resolution of the identity for the rotated probe state: n q n = ψ θ | n |n n| |ψ θ = 1 1 . (5.18) The classical Fisher information from the projections q n onto Bloch coherent states takes the form where we have averaged over Bloch coherent states. In this limit, measuring a sufficient number of projections lets one approach the QCRB, regardless of probe state, and a necessary condition for saturating the QCRB is the vanishing of the imaginary term in (5.19).

Metrological power of Majorana constellations
The geometrical structure of a quantum state informs its usefulness for rotation sensing. This is because rotations acting on a quantum state rigidly rotate the latter's Majorana constellation, making the rotation sensing problem equivalent to the distinguishing between constellations' orientations. We can see geometrically why some states are not useful for general rotation sensing, some states are useful for estimating rotations about known axes as in phase estimation, and others are useful for estimating rotations about unknown axes. We begin with SU(2) coherent states. Their constellations only possess two angular degrees of freedom, because they have a continuous rotational symmetry about one axis, so can never be used for sensing all three parameters of a rotation. This is reflected by the properties of the sensitivity covariance matrix (which we write for the particular case of a coherent state located at the north pole) The matrix has determinant zero, its condition number is formally infinite, and the trace of its inverse diverges. Next, we turn to eigenstates of angular momentum projection operators. Again, the associated constellations only have two angular degrees of freedom, which is why it cannot be used to simultaneously estimate all three rotation parameters. The sensitivity covariance matrix seems more useful than C(|JJ ), but it still cannot be inverted. It is well known that NOON states are useful for estimating rotations about known axes. Their constellations have an extra angular degree of freedom that yields a sensitivity covariance matrix Now we see that the matrix has a nonzero determinant J 4 /4, its condition number is J/2, and the trace of its inverse is (4J + 1) /J 2 . These states are useful for estimating rotations around a particular axis, evidenced by their geometrical representation and the C 33 component, but their performance for arbitrary rotations gets poorer with increasing J (as seen by the increasing condition number) and leads to a total uncertainty that only scales as O (1/J).
Bloch cat states also have a diagonal sensitivity covariance matrix, just like for NOON states. Still, the determinant decreases monotonically and the trace of the inverse increases monotonically with |Θ − π/2|, so we see that these Bloch cat states monotonically span the transition between the performance of NOON states and the inferior angular momentum eigenstates.
Majorana constellations that extend beyond a single great circle to make use of the three-dimensional structure of the unit sphere are instrumental to optimizing estimation of rotations about arbitrary axes and of all three rotation parameters. We can use this geometrical understanding to concoct other useful states for metrology. For example, to improve the usefulness of NOON states we can consider supplementing the stars on the equator of the unit sphere with some stars at the north and south poles: These now have 2m points spread about the equator of the unit sphere, supplemented by J − m points at each of the north and south poles (we take m > 1/2 without loss of generality to ensure that the states have vanishing angular momentum projections). How can we determine the optimal number of stars to put at the poles and along the equator? We know that these states will again be useful for estimating rotations about arbitrary axes, with most of their usefulness being for rotations about a single axis, when m ∼ J. We can determine the optimal balance through the sensitivity covariance matrix: We see that decreasing m from J while retaining m = O(J) makes the matrix more stable to inversion and the total uncertainty decrease, relative to NOON states. The optimal value of m 2 = 1 3 J(J + 1) yields a sensitivity covariance matrix C(|ψ balanced ) → 1 3 J(J + 1)1 1, the same as for the Kings; this simultaneously maximizes the determinant, minimizes the condition number, and minimizes the trace of the inverse. This optimal value of m can only be obtained for some values of J [2J ∈ (6,25,96,361,1350,5041,18816,70225,262086, 978121 · · ·)], making those states Kings, and for other values of J these states can only approach the Kings. This geometrical approach shows how rotation sensing performance can be improved by adding three-dimensional structure to the Majorana constellation.
The states with the best geometries are the Kings of Quantumness. These have Majorana constellations with discrete symmetries about multiple axes, whose highly distributed stars yield the sensitivity covariance matrix This has a determinant that scales with O J 6 , a perfect condition number of 1, and the trace of the inverse scales with O(1/J 2 ). From all of these perspectives it is clear that they have more metrological power for estimating arbitrary rotations than any of the other states considered in the context of phase estimation.
As an extension to these states, we can consider states that do not live in a single Hilbert space H J , as such states have been used in phase estimation, and compare their  (4,9,14,19) for two-mode photonic states with a coherent state in one mode and a squeezed state in the other. All of the constellations in each photon-number subspace lie on the same longitudinal great circle. The ratio of the coherent state strength to the squeezing parameter is α 2 /λ = 4 and a change in its phase merely rotates the constellation. The most sensitive states have |λ| ≈ 1; this is reflected in the constellations being the most evenly distributed for photon-number subspaces near |α 2 /λ|. performances to their geometrical representations. The Majorana representation must be extended to a set of nested spheres, one for each H J , and the information about relative phases between a state's components in each Hilbert space is lost. To make the connection to angular momentum, we represent each state by its occupation of two harmonic oscillator modes annihilated by a and b, and use the Schwinger map [90] The Greek index µ runs from 0 to 3, with σ 0 = 1 1 and {σ k } (k = 1, 2, 3) are the Pauli matrices. Note carefully that J = N/2, where N = a † a + b † b is the operator for the total number of excitations. For example, we can begin with the two-mode coherent states |α, β ∝ exp αa † + βb † |vac , (6.8) where |vac denotes the two-mode vacuum. In each subspace H J , this state is represented by 2J stars at the same angular location: the polar angle is always 2 tan −1 |β/α| and the azimuthal angle is always arg(β/α). These states have Similar to the Kings, the condition number is 1, making these states stable to matrix inversion. However, for average energies |α| 2 + |β| 2 and 2J, the two-mode coherent states scale much more poorly with energy in the context of determinant and total uncertainty. This can be interpreted from the geometry: the coherent states do not possess much angular information, because all of their stars lie on the same axis, but the hidden relative phase information can be used to determine rotations around arbitrary axis with some precision, because the states are not eigenstates of an angular momentum projection operator. The previous example is readily extended to other states useful for phase estimation. It can be shown that all of the stars in all of the subspaces lie along the same great circle for many such states, explaining why these are sensitive to estimating rotations about particular axes. For example, the two-mode squeezed states have half of their stars lying at opposite poles from each other in each subspace, giving more angular resolution than coherent states alone. A more intricate case is that of a state that is coherent in one mode and squeezed in the other: These states have Majorana polynomials ψ(z) proportional to 1 F 1 (−J; 1/2; −z/4) when J is an integer andz 1 F 1 (−J; 3/2; −z/4) when J is a half-odd integer, wherẽ z = −2α 2 z 2 /λ, α is the strength of the coherent state, λ = ξ |ξ| tanh |ξ|, and ξ is the strength of the squeezed state. Since these confluent hypergeometric functions all have J real rootsz k < 0, all of the Majorana stars lie about the same great circle, whose azimuthal angle varies with the phases of α and ξ. We plot the Majorana constellations in a number of subspaces for such a state in Figure 4, where it is clear that the points are all distributed about the same longitudinal line. The sensitivity covariance matrix is readily calculable; for simplicity we take α 2 ξ * to be real (and positive), making the matrix diagonal. If the average number of excitations in each mode is large, namely J a = |α| 2 /2 and J b = sinh 2 (|ξ|)/2, we find which is useful for estimating rotations about axes in a specific plane (determined by the relative phase of α and ξ) but not axes in all three dimensions. It has been wellestablished that the most sensitive such states have J a = J b 1 [91], for which |λ| ≈ 1; this is substantiated in the J 2 component of C(|ψ c+s ) and reflected in the constellations for subspaces with photon numbers |α| 2 being the most evenly distributed, per Figure 4. This spreading about a single great circle, which holds for all of the states in, for example, all of the states studied in Ref. [92], shows that states from phase estimation are not immediately tailored to multiparameter estimation.

Conclusions
Quantum sensing has become a distinct and rapidly growing area. The case of rotations is the epitome of how quantum properties can boost sensitivity and precision. Certain probe states are exceedingly sensitive to rotations, with the geometry of a probe state dictating its usefulness for rotation estimation. We have made an extensive use of Majorana's beautiful stellar representation: being purely geometric, it allows quantum mechanics to be reconciled with our physical intuition in the classical realm. This has allowed us to introduce in natural way a number of measurement schemes for simultaneously estimating all three parameters of a rotation and we expect these to be useful in years to come.
Singularities are relevant to QFI in at least two ways. First, the density matrix may be singular, as is the case for pure states, and its rank may vary with the parameters being estimated. When this happens, the QFI suffers a discontinuity that makes it differ pointwise from the Bures metric [54,55,93].
More physically, the QFI itself, and its matrix extension, is singular when a state cannot locally distinguish between some values of the parameters. This has been studied for classical FI from a Riemannian geometric perspective, where one can use the Moore-Penrose pseudoinverse [94] to estimate a subset of the parameters with a modified CRB [95] (although that may constitute an overly optimistic lower bound [96]). The quantum version should proceed along the same lines: one should first identify the set of parameters that are independent, then calculate the corresponding QFI matrix for those parameters [97,56]. In this Appendix we discuss physical scenarios that lead to singular QFI matrices and review the implications for metrology, using rotation estimation as an illustrative example.
The QFI matrix depends on the set of parameters through the evolved state and the semilogarithmic derivatives. For the QFI matrix to be singular, there must be a basis, and therefore a parametrization, in which an entire row and column vanishes. We can analyze this occurrence in the case of pure and mixed states.
For pure states, an entire row and column of the QFI matrix vanishes if and only if |∂ λ ψ ∝ |ψ for some parameter λ. For full rank density matrices, an entire row and column vanishes if and only if the weights p i of each eigenstate |i of the density matrix are independent of some parameter λ and all of the eigenstates satisfy the pure states' condition |∂ λ i ∝ |i . For sub-full rank density matrices, the latter condition is modified to only hold within the range of the density matrix: either the above condition ∂ λ |i ∝ |i holds or the rank of the density matrix is changing. Since the rank of the density matrix can only change while the weights of its eigenstates are changing, a singular QFI matrix implies the presence of a parameter for which each eigenstate has the dependence The condition on each of the eigenstates holds a) trivially, when |i is explicitly independent of λ, and b) when |i 's global phase depends on λ. The differential equation seems to imply that the local exponential growth or decay of an eigenstate with some parameter leads to a singular QFI matrix, but this is precluded due to the concomitant change in the weights of the eigenstates. We will see the recurring connection between singular QFI matrices and global phases throughout our examples.
A QFI matrix can be singular at a local set of coordinates for any state; alternatively, it can be singular for all sets of coordinates for certain states; and, finally, it may be singular for all sets of coordinates for all states. The first scenario may correspond to a coordinate singularity, the second typically heralds a shortcoming of the probe state, and the third often signifies an unseen interdependence between the parameters in question. We explore each of these scenarios in the context of rotation and related measurements.
Consider estimating the azimuthal component of a rotation axis when the polar coordinate is 0 or π; i.e., at the north or south pole. Similarly, consider estimating the difference between the first and third Euler angles in a ZY Z configuration when the second Euler angle vanishes. These points exemplify coordinate singularities, at which one of the three rotation parameters is undefined, and so no state can be used to estimate all three because ∂ λ |ψ = 0 for one of the parameters. Coordinate singularities can be identified by the QFI matrix being singular for all states and becoming invertible when any of the parameters is slightly perturbed or when a different parametrization with the same rank is considered. Within the original coordinate system, one can find the set of defined coordinates by the range of the QFI matrix; e.g., diagonalizing the matrix for an Euler angle parametrization shows that the null vector is associated with estimating the difference between the first and third angles [87]. Alternatively, one can find a new parametrization using a transformation whose Jacobian determinant squared cancels the singularity in the determinant. This is exemplified by switching from a spherical to a Cartesian coordinate system, for which the former's QFI matrix has a determinant proportional to θ 4 sin 2 Θ and the Jacobian governing the transformation has a determinant proportional to θ −2 csc Θ. Linear combinations of the new parameters will also all be estimable.
Conversely, when estimating the separation between two Gaussian sources, their relative intensities, and the centroid position, the QFI matrix diverges as the separation vanishes [98]. This singularity is present for all states at this local point, regardless of their coherence properties [99], implying a coordinate singularity. The rub is that a Jacobian erasing the singularity requires a strange coordinate system that is not conducive to measurement: for example, a transition from the separation coordinate to the square of that coordinate.
There are some states for which the QFI matrix is always singular regardless of parametrization. This is the case when estimating all three rotation parameters with classical states; Bloch coherent states |n cannot be used for estimating rotations about unknown axes. That Bloch coherent states are only sensitive to two rotation parameters is apparent from an Euler angle parametrization: by choosing the first rotation axis to coincide with the direction of the spin, the first Euler angle imparts a global phase on the state, making |∂ λ ψ = −i |ψ . One can again find this irrelevant coordinate by diagonalizing the QFI matrix for a classical state. It turns out that the only pure states with this issue are eigenstates of angular momentum operators J · n, as they are the only ones for which the covariance matrix Cov |ψ (J i , J j ) is singular, and this is the case regardless of parametrization. This immediately shows the advantage of using nonclassical states: they enable the simultaneous estimation of all three parameters of a rotation.
Finally, some parameters were never meant to be measured. For example, in the absence of a reference phase, one cannot detect an average phase shift among all of the branches of a linear optical network. Diagonalizing the QFI matrix shows that only relative phases are estimable, and this holds true regardless of the probe state and of the values of the parameters being estimated; as foreshadowed, the average phase imparts a global phase to the state [100]. Similarly, the QFI matrix is always singular when the parameters being estimated are not linearly independent. This scenario shows the necessity for new physics, such as the presence of an external reference phase, or new thinking, such as determining the true number of independent parameters, in order to estimate any of the parameters.
In all cases of singular QFI matrix we require a presciption for how to proceed with estimation. The trick detailed in Ref. [97] advises finding a parametrization in which the matrix is block diagonal and some maximal block is invertible, then proceeding to estimate the parameters corresponding to that block. The parameters can be inferred from a diagonalization of the QFI matrix, but we stress that a change in parametrization is not necessarily the same as a change of basis, as the corresponding Jacobian matrices need not be unitary. Moreover, one must avoid the use of finely tuned probe states for which preparation errors rid the QFI matrix of its block diagonal structure, as these cases require reintroduction of so-called nuisance parameters to the analysis. Singularities in the QFI matrix are essential to understanding the physical limitations of an estimation problem.