Light-ray operators, detectors and gravitational event shapes

Light-ray operators naturally arise from integrating Einstein equations at null infinity along the light-cone time. We associate light-ray operators to physical detectors on the celestial sphere and we provide explicit expressions in perturbation theory for their hard modes using the steepest descent technique. We then study their algebra in generic 4-dimensional QFTs of massless particles with integer spin, comparing with complexified Cordova-Shao algebra. For the case of gravity, the Bondi news squared term provides an extension of the ANEC operator at infinity to a shear-inclusive ANEC, which as a quantum operator gives the energy of all quanta of radiation in a particular direction on the sky. We finally provide a direct connection of the action of the shear-inclusive ANEC with detector event shapes and we study infrared-safe gravitational wave event shapes produced in the scattering of massive compact objects, computing the energy flux at infinity in the classical limit at leading order in the soft expansion.


Introduction
Event shapes or weighted cross-sections describe the distribution of outgoing particles and their properties such as charge or energy-momentum in a particular direction as determined by a detector. They are important tools for analyzing jets produced in collider experiments [1][2][3]. In the '90s, Korchemsky and Sterman 1 discovered that such event shapes can be equivalently described by the expectation value of some non-local operators -line integrals of local conserved currents -inserted at different locations on the celestial sphere [6,7]. While these objects are well-known in perturbative QCD, where they are used to study hadron scattering at high energies [8], it is only recently that Maldacena and Hofman initiated a systematic investigation of such operators in conformal field theory [9]. It is now well understood how to extract event shapes from OPE data and symmetries in conformal field theories (especially for N = 4 SYM) [10,11]. Moreover, the connection between energy event shapes and the ANEC operator establishes bounds on the a and c coefficients characterizing the conformal anomaly [9,12,13]. Of course, such developments led to increased interest in CFT light-ray operators which culminated in the systematic analysis of [14][15][16][17].
We can also study the correlation between event shapes by considering more than one detector (ie., multiple insertions of the corresponding non-local operators). Of particular interest is the so called energy-energy correlator, which is an infrared-finite observable that can be measured in collider experiments and used to provide a measure of the QCD coupling constant. The leading order QCD prediction for energy-energy correlators was computed in the late '70s [2,3]. More recently, analytical expressions at NLO [18,19] as well as numerical expressions at NNLO [20] have become available. In N = 4 SYM, there are even analytic results at NNLO [21].
The recent discovery of gravitational waves [22,23] started a new era in multi-messenger astronomy. A powerful technique to construct bank of template waveforms for the detectors is the effective one-body (EOB) map [24,25] developed in a seminal paper by A. Buonanno and T. Damour [26], which combines both analytical and numerical (when combined with numerical relativity (EOB-NR)) data for the description of the bound state problem in general relativity. It was pointed out by Damour [27,28] that the scattering problem and its post-Minkowskian expansion contains useful information which can help to construct the effective one-body hamiltonian for the inspiral phase. More recently the relation between bound and unbound orbits has been better clarified, at least in the conservative regime [29][30][31].
Scattering amplitudes and on-shell techniques have proven to be very useful in order to study the conservative regime for compact binaries. In an EFT approach [32] such objects are treated as point particles with possible additional finite size effect operators. Taking the classical limit from quantum scattering amplitudes requires to understand which type of coupling gives the classical spin dynamics [33][34][35], which was referred as "minimal coupling" in the literature [36,37]. With the matching with a suitable EFT potential [38] and combining several techniques like double copy and unitarity, it was possible to extract for the first time the conservative binary hamiltonian 2 at G 3 order for the spinless case [41,42] and at G 2 order for the spinning case [43]. This result has passed many non-trivial checks at the quantum [44] and at the classical level [45][46][47]. There has been some attempts to study radiation in this context both using classical worldline calculations [48][49][50][51][52][53][54] and double copy techniques scattering amplitudes in the classical limit [55,56]. In addition, several classical observables like the classical impulse [57], the radiation kernel [57] and the classical spin kick [35,58] were defined in full generality from from the classical limit of matrix elements of appropriate quantum operators.
Here we would like to understand how the energy carried away by gravitational waves, which is collected by a detector localized in a directionn, is described at the quantum level. We will see that the corresponding operator is related to the Isaacson effective stress tensor [59,60] by the integration over the retarded working time of the detector, and its action can be computed explicitly using standard techniques in the asymptotic expansion. Using the Bondi gauge framework for asymptotically flat spacetimes, we will show that the energy flow operator can be also written in terms of the Bondi news squared term and this provides an extension of ANEC operator at null infinity to a shear inclusive ANEC. With our detector interpretation, the latter naturally provides the sum of energies of all massless quanta (i.e. radiation) emitted in a directionn when acting on on-shell states. One can then generalize the standard notion of event shapes used in QCD to the gravitational case, which we will study in detail in the soft kinematic regime by taking the classical limit and summing explicitly over all graviton contributions. Indeed in the very low-energy regime, a single coherent state provides the exact quantum state of radiation and taking the limit of the large number of quanta reproduce known classical limit results for the energy flux [61][62][63][64]. A systematic study of the low-energy gravitational radiation produced in the scattering of compact objects up to sub-subleading order in the soft expansion has already been performed in the literature, both in the quantum and in the classical context [65][66][67][68][69][70][71][72]. We will also study the correlation between different gravitational energy event shapes, especially in relation to the expected classical soft factorization.
Light-ray operators constructed from the components of the stress tensor are naturally related to the global Poincaré charges by turning the line integral in the light-ray definition into an integral over all of space. We will restrict our analysis to the subset of light-ray operators obtained by integrating components of the stress tensor over a null-sheet. These light-ray operators were recently studied in the context of unitary CFTs by Cordova and Shao [73] and shown, via symmetry arguments, to form a closed algebra. While it is not immediately clear that the same light-ray operators will form a closed algebra in a generic QFT, the underlying universality of the stress tensor algebra [74] hints at this possibility. The stress tensor algebra was also reconsidered more recently [75,76] to derive a "universal" effective light-cone algebra for CFTs in d > 2. Moreover, there is an intriguing relation with the BMS algebra [77,78] -and its extended version [79][80][81][82][83] -with light-ray operators [73] which is still worth exploring for light-ray placed at null infinity. Motivated by that, we will define a new system of light-ray operators for linearized gravity, which will naturally combine with the standard stress tensor definition in scalar and gauge theories providing a unified treatment of all massless particles at I ± . We will then compute in our approach the associated algebra of light-ray operators in 3 + 1-dimensional QFTs for massless particles of integer spin and compare with [73].
Preliminaries Let's consider the process of electron-positron annihilation in QCD: for a generic outcome of the scattering process e + e − → X, (1.1) we would like to understand the properties of the emitted particles as captured by one (or many) physical detectors located at spatial infinity. To make contact with the structure at null infinity, it is convenient to work with flat null coordinates for Minkowski spacetime so that the null boundaries are located at r → +∞ while keeping (u, z,z) fixed. On this hypersurface (z,z) are stereographic coordinates on the celestial sphere. This set of coordinates correspond to wrapping up the transverse spatial coordinates of the light-sheet at infinity via the transformation Instead of considering the total cross section for the process (1.1) where p 1 and p 2 are the incoming momenta of e + and e − and k X is the sum of the outgoing momenta, we introduce a less inclusive observable by defining a weight w(X) such that which is the so-called "weighted cross section". Different choices of the weight factors w(X) can tell us about different properties of the distribution of outgoing particles, like the flow of quantum number or of energy and momentum in a particular direction on the celestial sphere.
Of particular interest is the to the choice of w(X) corresponding to the energy flow. Given an external on-shell state |X = |k 1 . . . k n of n massless particles with k µ where k µ i = (E i , k i ) is the 4-momentum of the individual particles and Ωk i = k i |k i | is the solid angle in the direction of k i [10,11]. The weighted cross section measures the distribution of energy in the final state that flows in the direction ofn, and as such we can interpret w E (X) as an eigenvalue of the average null energy operator (ANEC) at spatial infinitỹ where the energy momentum tensor is always understood to be normal ordered and t is the physical working time of the detector. Looking at the energy flow operator (1.8) it is easy to generalize it to a momentum flow operator (which is well-known in the jet literature [84,85])P whose action on a state gives the linear momentum flowing the particular directionn One can define also the expectation value of the ANEC operator at infinity, namely the 1-point energy event shape where it is always understood thatẼ(n) acts on the set of outgoing states |X inserted via the completeness relation, as usual in the Schwinger-Keldysh formalism. The radiation due to the scattering of two massive objects is captured by two detectors located at spatial infinity in the directionsn 1 (z 1 ,z 1 ) andn 2 (z 2 ,z 2 ).
In the case of two (or more) detectors, as previously mentioned, we can also define the energy-energy correlator which is related to the 2-pt Wightman function of two ANEC operators at infinity (1.8) inserted at different points on the celestial sphere It is quite interesting to analyze this 2-pt correlator for the gravitational radiation in the classical limit, where in the most general case a superposition of coherent states represents the radiation at the quantum level.
In flat null coordinates the momentum flow operator can be written directly in terms of the energy flow where u is the retarded time and identifyn with the related coordinates on the celestial spheren ⇔ (zn,zn). (1.14) Here, we emphasize that we use E(n) in place ofẼ(n) for technical convenience. One can work in Bondi gauge where the detector would register the energy flow as given byẼ(n) (i.e. the standard local energy factor). We will do this later on for gravitational event shapes.

Light-ray operators in QFT
In this section, we define the family of light-ray operators under consideration and make contact with their physical significance and in particular with the corresponding surface charges of null-sheet symmetry generators. Light-ray operators depend implicitly on a choice of null-sheet and we choose future null infinity. In flat-null coordinates, our family of light-ray operators is defined by We will collectively denote these operators by L(n) = {E(n), K(n), N z (n), Nz(n)}. This family of light-ray operators can be understood as surface densities 3 for the corresponding conserved charges of the future null boundary of Minkowski spacetime. To see this, consider an affine transformation of the null surface generator [86] δu = A + B u, (2.2) where A and B are arbitrary constants. Since lim r→∞ r 2 T uu (u, r, z,z) corresponds to the generator of null diffeomorphisms, we can integrate it to obtain either the generator of null which is an equivalent definition of d 2 z K(z,z). The latter is called "boost mass" M boost in the literature and it is usually defined for null-like horizons 4 using the boost Killing vector of a causal wedge in Minkowski spacetime [86][87][88]. In a CFT, this provides a direct connection to the concept of modular hamiltonian for null-like horizons [89]. While for empty Minkowski spacetime the boost mass is exactly zero, for more general non-asymptotically flat spacetimes it might be relevant for a proper formulation of the first law of black hole thermodynamics [90]. In particular, it turns out that while the ADM mass M ADM measures the total "monopole" distribution of the matter in the spacetime, the boost mass M boost measures the total "dipole" moment of the mass distribution at infinity [90].
Since we are only considering QFTs on a flat Minkowski background, we will not delve deeper into these topics. However, we will find the monopole versus dipole distinction relevant for understanding the action of K(n) on on-shell states. From now on, we will call K(n) the boost energy (surface) density. In the same spirit, the ANEC at infinity can be called null energy (surface) density.
Following a similar line of reasoning one finds that N z (n) and Nz(n) are associated to the local vector fields ∂ z and ∂z, i.e. with the angular momentum flux charge Therefore N z (n) and Nz(n) can be thought as angular momentum flux (surface) densities. More generally, this family of non-local light-ray operators appear naturally in Einstein equations solved near null infinity and in particular in (BMS and Poincaré) balance flux laws [91], where they represent the radiative contribution to the fluxes. 4 More generically, choosing a spatial slice of the null surface u 0 = Λ(z,z), a canonical "boost energy" can be defined [86] du (u − Λ(z,z)) lim r→∞ r 2 Tuu(u, r, z,z), (2.5) which turns out to be connected with the notion of "area operator", providing thus a link to the concept of generalized entropy and to the generalized second law.
Here we interpret physically the insertion of our light-ray operators as an insertion of a physical detector on the boundary of Minkowski spacetime. We will be interested in how these operators act on on-shell states and we will then compute their expectation values in perturbation theory, as in the standard event shapes literature.

The Spin-0 light-ray operators
Here, we consider a theory of self-interacting massless scalars with a potential Q(φ) that is polynomial in the fields and without derivative interactions. We derive explicit expressions for our family of light-ray operators in terms of the scalar creation and annihilation operators.
The Hilbert stress tensor T matter for our theory of scalars is simply The scalar stress tensor components which are relevant for the light-ray operators are 5 The next step is to substitute the leading term in the large r expansion of the scalar field into the expressions above. However, to do this, we need to impose some boundary conditions on the field. Since we require the energy flux across I + to be finite, we impose the following falloff condition for the scalar field 6 φ(u, r, z,z) = φ (0) (u, z,z) r Far away from interactions we can safely use the free field mode expansion (3.4) and evaluate the field with the large r saddle-point estimate. Physically, r is the distance from the detector to the region where the interaction is localized and it should be much larger than the inverse of the smallest momentum scale that appears in the S-matrix [85]. Note that in our coordinate conventions the canonical commutators read One can easily work in a more general setting with derivative interactions. However, such terms will drop off due to their scaling in 1 r in the limit that the null-sheet is pushed out to null-infinity. Instead, we prefer the interaction terms drop off directly from the structure of stress tensor components since this more closely resembles scalar light-ray operators in the bulk (see appendix C). 6 Soft scalar modes are beyond the scope of this work, so we do not allow the scalar field asymptotic expansion to have a constant piece around I + ± .
The large r limit of such a Fourier integral is controlled by the exponents of the exponentials where we use as our parametrization for the coordinates and on-shell momentum [92,93]. The only saddle point is at (z,z) = (w,w). Making the following change of variables (y 1 , y 2 ) = 1 2 (z +z), − i 2 (z −z) brings the exponentials into Gaussian form, which integrates to Substituting into the mode expansion for φ(x) yields the leading term in the large r limit When derivatives in z (orz) acts on field operators expressions, one can write 7 for the saddle point estimate of ∂ z φ (0) where the derivative acts on the creation and annihilation operators localized around (w,w) = (zn,zn).
Inserting the saddle point estimate of φ(x) into the definition for the scalar ANEC operator at infinity, we find where we have set to zero all contributions proportional to δ(ω p 1 + ω p 2 ) since the only cases that these constraints can be satisfied correspond to ω p 1 → 0 ± and ω p 2 → 0 ∓ , which are unphysical. The action of E scalar (n) on an on-shell state |X = |p 1 ...p n gives which is natural for an observer located at spatial infinity who is measuring the energy flux along the retarded time u in flat null coordinates. Note that on-shell scattering states are eigenstates of the ANEC operator in the detector limit whose eigenvalues are weight functions similar to equation (1.7).
Inserting (3.9) into the definition of the boost energy density flux K scalar (n) yields It is convenient to change variables where ω + ∈ [max{−ω − , ω − }, +∞[ and ω − ∈] − ∞, +∞[ is a convenient slicing of the integration region in the new variables. 8 In the new variables, 9 where we have used the following distributional identity (for n = 1) The derivative acting on the energy delta function might seem troubling at first sight.
To develop more intuition, we study the action of this operator on on-shell states |X = |p 1 ...p n : Then one can solve the ω + integral where we have relabelled ω − as ω. It is also enlightening to compare the simplest matrix elements of E scalar (n) and K scalar (n) The original frequencies are given by ωp 1 = (ω+ +ω−) and ωp 2 = (ω+ −ω−) and the integration measure becomes dωp 1 dωp 2 = 2 dω − dω + . 9 Here we have assumed, as before, that the other contributions dropped because they are unphysical.
It is well-known that in QFT the single contraction q|p must be interpreted in a distributional sense; in order to properly define such objects one needs to smear them with a well-behaved function (see [94,95] and references therein). In the S-matrix context, one usually considers outgoing states of the form 10 where φ p (p) is a suitable real momentum wavefunction localized around p. Note that Looking at (3.21), the interpretation of K scalar (n) is also much more clear: while δ(ωq − ωp) can be considered the limit of a suitable localized function around ωq = ωp, it turns out that is probing a dipole-like pattern around ωq = ωp (see 2a and 2b respectively).
(b) Dipole profile centered around ωq = ωp, which converges to δ (1) (ωq − ωp) in the limit of zero width. Figure 2: Comparison of the expectation values of ψ q,out | E scalar (n) |ψ p,out and ψ q,out | K scalar (n) |ψ p,out as a function of the energy ω.
The final operator we consider is N z,scalar (n), which after using (3.9) and (3.10) is More compactly, is the standard orbital angular momentum operator. A similar analysis can be repeated for the antiholomorphic component by flipping z ↔z.

The Spin-1 light-ray operators
Having found explicit expressions for the light-ray operators for massless scalar theories, we now derive similar results for the light-ray operators for massless spin-1 theories. Using the Hilbert stress tensor T matter for pure gauge theories, one obtains Here, the SU (N ) generators are chosen to be hermitian [T a , T b ] = if abc T c and are normalized according to Tr(T a T b ) = 1 2 δ ab . To discuss the large r limit of the stress-tensor, we need the large r limit of the gauge field, which requires making a gauge choice, as well as, specifying boundary conditions. We will adopt the radiation gauge [96] A u = 0, and we require the gauge field components to satisfy the following falloff conditions at infinity which are needed to have a finite energy and angular momentum flux at I + . The stress tensor components relevant for our family of light-ray operators (4.2) are up to terms of order O(1/r 3 ).
We assume that we can work perturbatively with energies above the confinement scale Λ > 0 and can formally talk about asymptotic gluon states. We use the following parametrization and projection of the polarization vectors and the canonical commutation relations Starting with the gluon-ANEC operator at null infinity, we find where the saddle point estimates for (4.12) Analogous to the scalar case, on-shell gluon states |X = |{p 1 , σ 1 }, ..., {p n , σ n } are eigenstates of the gluon-ANEC at infinity The gluon boost energy density operator K gluon (n) follows straightforwardly (4.14) These operators are especially simple; the only fundamental difference between {E gluon , K gluon } and {E scalar , K scalar } is the sum over helicities.
On the other hand, the spin-1 angular momentum flux operators are complicated by boundary terms, which come from integrating by parts in u. Explicitly, where we have used the leading u-equation of motion It is worth remarking that the non-abelian contribution appears through the boundary term (4.15); this is similar to the case when the light-sheet is placed in the bulk (see appendix C). For the first contribution to the angular momentum flux density, we use the identity just to isolate the contribution coming from the total derivative. Regarding the saddle point This gives the following representation for the angular momentum flux density • An orbital angular momentum contribution This piece is expected from the spin structure of the point particle stress tensor [97]. Indeed the covariant spin matrix reads [98] (with σ helicity) and the component S zz = 1 r 2 (∂ z x µ )(∂zx ν )S µν = i σ 2 enters in the total hard angular momentum operator as shown in appendix A of [99].
Please note that the splitting of angular momentum into an orbital and spin part is not gauge invariant, and it is done here only for convenience in analyzing different terms. The boundary terms in (4.15) are associated with soft contributions 11 . Using the equations of motion, we can write This makes it clear that for pure Maxwell theory (which is conformal) this term vanishes identically. Using the saddle point we obtain where in the last line we recognize the gluon number density operator which contributes to the hard part of the non-abelian (large gauge) charge [93,100]. Defining the boundary value of the field [93, 96] the total boundary contribution to the angular momentum flux becomes This contribution is soft, since it contains the soft mode C a z (z,z) 12 . It is also contains the hard gluon charge, which might seem strange at first. We can gain some insight by turning on a matter source j µ , which modifies the u-equation of motion and therefore the final density of angular momentum flux du lim r→+∞ where is the density of the hard color charge. These terms were first found by Ashtekar et al. [101,102] in the abelian spin 1 case (coupled with a matter distribution) and they are related to mixing of the so-called radiative modes (i.e. our A ) and of coulombic modes (i.e. the hard charge distribution, potential modes). A similar expression holds for Nz ,gluon (n) which can be obtained by exchanging z ↔z and flipping all helicity terms (+) ↔ (−).

The Spin-2 light-ray operators
It is well-known that the equivalence principle prevents a universal definition of a stress tensor for the gravitational field. One can try to define various nontensorial objects called "pseudotensors" which when linearized transform as Lorentz tensor. However, all such definitions are not gauge invariant (a notable example is the Landau-Lifshitz pseudotensor [103]). In this paper, we will only consider linearized gravity, i.e. we split the metric tensor as g µν = η µν + h µν and assume |h µν (x)| 1. In the physically relevant situation of gravitational waves propagating in an asymptotically flat region, a suitable notion of a gauge-invariant stress tensor for long-wavelength modes can be defined using the so-called Isaacson averaging procedure [59,60]. In our conventions, the effective stress tensor for gravitational waves is [104] T eff,GW where · stands for the average over short-wavelength graviton modes. 13 It is possible to extend the Isaacson construction of an effective gravitational wave stress tensor to general classical theories of gravity [105]: for example, the leading order contribution at infinity is proportional to the Ricci mode for f (R) theories [106], conversely the spectrum is unchanged for most higher derivative theories of gravity since higher derivative terms usually drop off [105]. It is worth stressing that this stress-tensor is defined canonically a la' Noether, and therefore it represents only the physical energy and the momentum flux gravitational wave contribution whereas the angular momentum requires a different analysis.
Asymptotically flat geometries have the following large r asymptotic expansion near I + in the Bondi gauge [77,78] where the leading contribution comes from Minkowski spacetime and subleading contributions encode deviations from the flat background. Bondi coordinates are defined as is the Bondi angular momentum aspect. Here we use the conventions of Pasterski-Strominger-Zhiboedev [99], which differs from Barnich-Troessaert [81], Hawking-Perry-Strominger [107] or Flanagan-Nichols [108]. Since different conventions change the definition of the light-ray operators, this is of particular relevance here. 14 Spaces that admit expansions of the form (5.2) are called Christodoulou-Klainerman (CK) spaces in the literature [110,111]. Physically, the shear tensor encodes gravitational radiation analogously to Maxwell field strength 15 In these coordinates, where κ = √ 32πG [112]. The time-time component of the Einstein equations (5.2) gives the Bondi mass loss formula For CK spaces, the asymptotics of the Weyl curvature component Ψ 0 2 implies the following conditions [110] m B (v, ζ,ζ) 14 See also [109] where some of the relations between different conventions are spelled out in detail. 15 In light-cone gauge, Av = 0.
where C(ζ,ζ) labels inequivalent BMS vacua related to each other under action of the supertranslation mode. After integrating over v from −∞ to +∞, we get [112,113] (ζ,ζ) = 0 for flat Minkowski spacetime. Using supertranslation invariance, one can choose a Bondi frame [110] where C ζζ (−∞, ζ,ζ) = 0 but in general the final vacuum will be non-trivial C ζζ (+∞, ζ,ζ) = 0. Focusing on Bondi news squared term we see that it agrees with the (light-cone) time-time component of Isaacson effective gravitational wave stress tensor at infinity T eff,GW vv in these coordinates (albeit this is defined with the implicit averaging over short wavelength modes). The action of (5.8) on an on-shell n-graviton state |X = |{p 1 , σ 1 }...{p n , σ n } is where we used the saddle-point estimate for the free field mode expansion for the linearized graviton mode The weight factor of the Bondi news squared term is given by the standard local factor E p in retarded Bondi coordinates, which is expected because here v = t − r is the standard light-cone time.
We can now come back to the equation ( where we have defined a suitable extension of the ANEC operator to a shear-inclusive ANEC at infinity which includes both matter and gravity contributions whose action on an on-shell n-particle state of different massless species α∈species |X α is (5.14) As we will prove in the next chapter on gravitational event shapes, expectation values of this operator will be positive-definite for a scattering process where we treat (perturbative) gravity as an EFT, similarly to what happens in the massless spin 0 and spin 1 case. Originally, the shear-inclusive ANEC operator was defined for complete achronal null geodesics (also called null lines) to cure the violations of the averaged null energy condition appearing for linearized graviton perturbations, see the appendix D for more details.
At the quantum level we can promote both sides of (5.12) to operators and act on an on-shell state composed of massless particles reaching null infinity (i.e. radiation): a direct calculation shows that where ∆M B = 0 in flat Minkowski picture at v → ±∞. This is the leading memory effect due to the radiation flux generated in the scattering process [114], which can be directly related to Christodoulou memory effect [115]. Our results are consistent with the memory being given by the (transverse traceless part of) soft factor [113,114] M ζζ = 1 once the normalization factor 1/(4πr) is taken into account. The angular momentum flux, requires the Einstein equations of the subleading terms in the Bondi gauge expansion [99] ∂ Like the shear-inclusive ANEC at infinity, all angular momentum flux contributions that are quadratic in the fields combine together into a single expression for all massless particles where the gravity contributions are expressed only in terms of the radiative data. The saddle point estimate of the first two terms of such contributions above gives the density of orbital angular momentum flux The remaining term gives the spin contributioñ Other conventions like Hawking-Perry-Strominger (HPS) differ from (5.18) in total derivative terms in ζ orζ of C ζζ N ζζ and CζζNζζ, which look like contact/spin terms under a saddle-point analysis similar to the one discuss before. While we are not attempting any kind of rigorous analysis of all these terms in different conventions here, our expectation is that the PSZ convention is the most natural 16 choice which gives the expected spin term of the graviton. This was explored further in [91], where the most general expression of the total angular momentum flux was analyzed and it was shown that (5.18) corresponds to the simplest case α = β = 0 in their notation.
To make contact with the system of coordinates used in the rest of the paper, we re-write the corresponding expressions in flat null coordinates 17 whose saddle point estimate gives : a † σ (ω pn )a σ (ω pn ) :, Note that the spin operator above differs from the spin-1 counterpart by a factor of 2 since the graviton has helicities σ = ±2. 16 Probably also because it makes direct contact with the Weyl tensor component Ψ 0 1 . 17 These expressions can be derived in other ways. For example, EGR(n) and KGR(n) can be deduced from the saddle point estimate of the Isaacson effective stress tensor in flat null coordinates.

On-shell light-ray algebra
In general, computing the commutator algebra of composite (non-local 18 ) operators is quite subtle unless completely fixed by symmetries. For example, with the Schwinger action principle [116][117][118] we can compute the structure of all commutators of the (covariantly conserved) stress tensor but the commutators of stress tensor components with spatialspatial indices; such commutators are model-dependent [74]. Since our family of light-ray operators contain either the time-time components or mixed time-spatial components of the stress tensor, it is conceivable that our light-ray operators satisfy a universal algebra [73,75,76] [L 1 (n 1 ), L 2 (n 2 )] = L 3 C L 1 L 2 L 3 (n 1 ,n 2 ) L 3 (n 2 ). (6.1) To make contact with physical detectors at infinity, we study the light ray algebra by taking the difference between the expectation value of Wightman two-point functions: Without loss of generality, we can set the |ψ above to be a single particle state. In the following calculations we will assume that canonical commutation relations are valid, motivated by the fact we will compute the algebra on a flat section of null infinity thanks to the choice of flat null coordinates.

Spin 0
The simplest commutator is P µ scalar (n 1 ) and E scalar (n 2 ): q| [P µ scalar (n 1 ), E scalar (n 2 )] |p = 0. (6. 3) The vanishing of this commutator -at all order in perturbation theory-is physically interpreted as the statement that measurements of energy and momentum in two localized directionsn 1 andn 2 are compatible at the quantum mechanical level.
The commutator [K scalar (n 1 ), E scalar (n 2 )] is slightly more non-trivial. Using q| [K scalar (n 1 ), E scalar (n 2 )] |p = δ 2 (zn 2 − zp)ω p − δ 2 (zn 2 − zq)ω q q| K scalar (n 1 ) |p , (6.4) 18 The fact that our operators are non-local does not guarantee the convergence of the OPE of such operators in the first place. This problem can be solved in CFT [15], but it is not clear how to generalize it to a general quantum field theory.

Spin 1
The spin 1 abelian case (i.e. Maxwell theory) does not differ much from the spin 0 case except for the presence of the helicity term. The additional helicity term commutes with : a † σ (ω pn )a σ (ω pn ) :, because it is proportional to the number operator which counts the difference between the number of helicity plus and helicity minus photons The commutation relation of N orb z with N orb z is exactly the same as the scalar case This is obvious from the fact that the N orb z,photon is simply the sum of two copies of N orb z,scalar N orb z,photon (n) = i d 3 p δ 2 (zn − zp) σ=±1 : a † σ (ω pn ) ↔ ∂ zn a σ (ω pn ) : . (6.14) Moreover it is clear that N spin z,photon (n 1 ) commutes with N spin z,photon (n 2 ) [N spin z,photon (n 1 ), N spin z,photon (n 2 ) = 0. (6.15) The remaining orbital-spin commutators take the form [N orb z,photon (n 1 ), N spin z,photon (n 2 ) = +2i∂ z δ 2 (zn 1 − zn 2 )N spin z,photon (n 2 ), [N spin z,photon (n 1 ), N orb z,photon (n 2 ) = −2i∂zδ 2 (zn 1 − zn 2 )N spin z,photon (n 2 ) + 2iδ 2 (zn 1 − zn 2 )∂zN spin z,photon (n 2 ). (6.16) The above commutators are quite a bit more subtle than the others encountered so far and we do find an ambiguity in their derivation. 20 However, any problems with the orbital-spin commutator disappear in the integrated charges since N spin z/z are total derivatives of z orz. 20 Since the N z/z do not contain integrals over the transverse directions, one must smear them with a test function. A more rigorous analysis would involve studying Wightman functions in position space with a careful prescription for sending the operator positions to the same light-sheet (see [119], for example) in addition to smearing with test functions.
If we turn on matter contributions, we find an additional (hard) interaction term which is mixing Coulombic and radiative data du A (0),a z j a,(2) u , (6.17) which clearly breaks the algebra because it does not conserve the particle number.
Regarding the soft non-abelian contribution which comes from the boundary term, we stress that as in the spin 0 case we are only interested in the hard contributions and therefore ignore the soft terms in our analysis, which requires a more delicate study of the soft non-abelian sector and its quantization as done by He and Mitra in [93]. We expect such terms to be related with the gluon soft theorem [96], and therefore to be relevant for the light-ray algebra of pure Yang-Mills in the soft sector: the non-commutativity of soft gluon limits for different helicities might give a non-trivial extension of Cordova-Shao algebra. We leave this very interesting problem for a future study.

Spin 2
Here similar comments to the ones made in the spin 1 abelian case apply, except for the fact that we have some freedom in the choice of the Bondi angular momentum aspect and therefore of the density of the angular momentum flux as discussed before. We find that with our definitions the algebra of light-ray operators for the graviton case is consistent with the complexified Cordova-Shao algebra. In particular, there are no contributions which are mixing coulombic and radiative data as expected [120]. Moreover, the convention adopted by Pasterski-Strominger-Zhiboedev [99] seems to be the most natural for the system of light-ray operators since it provides the standard helicity term which we would expect from a spin 2 point particle stress tensor.

Comparison with complexified Cordova-Shao algebra
It is straightforward to check that the commutators [E(n 1 ), E(n 2 )] = 0, [K(n 1 ), K(n 2 )] = 0, agree exactly with [73] since T uu = T −− 21 . However, in order to compare the other commutators, we need to compare how the components of the stress tensor at null infinity in our cooordinates differ from those in the standard light-cone components It is worth remembering here that δ 2 (Ωn 1 − Ωn 2 ) = 2δ 2 (zn 1 − zn 2 ).
Therefore, we see that a full detailed comparison with the (complexified) Cordova Shao algebra for the angular momentum flux requires understanding the following operator where there is an explicit mixing with the original contribution of the complexified version of N A at infinity. In general we have (r = x + ) lim r→+∞ r 2 T uu = lim match the complexified Cordova-Shao algebra. The only commutation relation which is new and differ from their result is [N z (n 1 ), Nz(n 2 )] 22 , due to mixing in (6.22).

Gravitational energy event shapes and S-matrix
In analogy to the standard event shapes considered in QCD jet physics, we explore gravitational event shapes, which are of current interest to the gravitational waves community. Our light-ray operators probe the underlying structure of gravitational radiation as registered by a detector placed at null infinity in the directionn. Of particular interest is the gravitational radiation emitted by the scattering of compact objects. We consider gravity as an effective field theory valid below the Planck scale [121,122] and work in the weak field regime where we can rely on perturbation theory. We assume, for simplicity, that such compact objects are spinless and point-like; that is, they can be described massive scalar point particles of mass m 1 and m 2 minimally coupled to gravity. The relevant amplitude for such processes is of the form where |{k j , σ j } j=1,..,M are the graviton momentum eigenstates and |p i are the momentum eigenstates of the two massive point particles.
In general, every radiation state in QFT can be expressed in terms of coherent states |α which provide an overcomplete basis [123,124]. Moreover in the case where the classical limit is well-defined the so-called Glauber-Sudarshan P-representation [123,125] of the radiation density matrixρ admits only positive weights P (α) ≥ 0. In particular in suitable circumstances for soft radiation [70,126] or in the eikonal limit [71,127,128] a single coherent state can represent the classical radiation at the quantum level. Therefore, as a first approximation, we focus on the soft kinematic region of the emitted gravitons in order to take advantage of the simple form of graviton soft factors and of the fact the quantum state representing the radiation is known. To be more precise, the detectors are assumed to have a lower energy resolution λ and the emitted gravitons to have total energy E wave with where Λ is a suitable UV cutoff which can be thought as the total energy of the process [64,70]. The classical limit requires the number of soft gravitons emitted to be large M → +∞, as dictated by the coherent state structure [65]: this implies that the limit is quite delicate in general, and it implies an infinite resummation over all graviton contributions. Please note that we use Bondi retarded coordinates instead of flat-null coordinates in this section, in order to express everything in terms of the standard energy weight factor.

The 1-point gravitational energy event shape
The 1-point gravitational energy event shape is simply the on-shell expectation value of the gravitational ANEC at infinity out|Ẽ(n) |out = in| S †Ẽ (n)S |in . (7.4) Inserting the completeness relation of on-shell states X |X X| yields which is always positive definite. The graviton emission amplitude for M soft gravitons of momenta k j and helicity σ j is then 23 At leading order there is no distinction between consecutive double soft limits and simultaneous double soft limits, so we don't need to assume any hierarchy in the energy of the soft gravitons.
We will be interested in the classical limit of this quantity, and the proper way to do so would be to use the in-in KMOC formalism [57]. In such approach, the limit is taken with a careful analysis of the wavepackets of the external massive fields which localize the particles on their classical trajectory as → 0. Here for simplicity we will work in the approximation ω b −1 , where ω is the typical graviton frequency of the emitted wave and b is the impact parameter of the two incoming particles: this will allow us to get a simple universal b-independent result which will be useful to analyze the properties of gravitational energy event shapes. In order to make transparent the power counting for the amplitudes involved, we will consider the tree level classical contribution for A GR (p 1 , p 2 ; p 3 , p 4 ) so that |A GR (p 1 , p 2 ; p 3 , p 4 )| 2 will always be of order G 2 in this section.
The leading contribution will be given by the single graviton emission amplitude, In the soft regime for the graviton we can isolate the graviton phase space integration, To make further progress we need to use the polarization sum identity for the graviton polarization vectors which reads σ=±2 σ µν (k) σ, * αβ (k) = 1 2 [P µα P νβ + P µβ P να − P µν P αβ ] , where q µ is a suitable reference momentum. The gauge dependent terms will cancel from our calculation using momentum conservation [64], so we do not need to worry about spurious contributions. Since the amplitude factorizes, we can actually perform the graviton phase space integration in a universal form Working in the center of mass frame 24 for the two massive particles we have which is infrared finite at the leading order G 3 . It is worth noticing that (7.11) is only the leading order contribution, and we still need to sum over all possible graviton emissions. Before looking at the classical limit, it is worth looking at how infrared divergences will 24 Here s = (p1 + p2) 2 is the center of mass energy.
cancel for the gravitational energy event shape in general. At order G 4 we have, The first contribution includes loop contributions related to virtual gravitons where we have regularized the integral with the ultraviolet cutoff Λ [64]. The second term instead has a divergent contribution which comes from the emission of one real graviton, since integrating over the 2-gravitons phase space gives where we have used the fact that which is valid as long as we can neglect the single graviton energies (suppressed by ) compared to the energy of the wave E k 1 , E k 2 E wave , as noticed in [70].
One can then perform the regularized sum in (7.12) which is infrared finite as expected. One can then repeat the same argument for the G 3+n -th order contribution It is worth noticing that this situation is in contrast with QCD case where the gluon energy event shape is infrared divergent and one has to sum also over quarks energy event shapes to get a well-defined (i.e. an infrared finite) answer [1]. The cancellation of infrared divergences in this setting is similar to what was found for the graviton emission cross section by Donoghue [133].

The 2-points gravitational energy event shape
The main object of interest here will be out|Ẽ(n 1 )Ẽ (n 2 ) |out = in| S †Ẽ (n 1 )Ẽ (n 2 )S |in . (7.18) By inserting a completeness relation of on-shell final states X |X X| between the two gravitational ANEC operators at infinity we get In details Contrary to the 1-pt energy event shape, here we get also contact terms which correspond to the case where the two detectors are aligned along the same direction. These terms are usually included in QCD energy event shapes and they are required in order to remove collinear divergences [2], but classical gravity does not have such type of collinear divergences [134]. We will thus restrict our attention to the generic casen 1 =n 2 : At this point we can repeat the integration as we did before to get As we can see, despite being a natural quantum observable the 2-point gravitational energy event shape factorizes, i.e. there is no degree of correlation between the two emissions. This is due to the soft limit we are considering and ultimately related to the coherent state structure of the radiation [71]. It would be very interesting to compute the non-trivial (connected) 2-point energy event shape in a quantum theory of gravity from the six point amplitude with 4 massive scalars and 2 gravitons, which should provide the leading infrared finite contribution to this observable.

Gravitational wave energy event shapes and classical factorization
So far we have discussed the quantum picture, where energy event shape correlators are a non-trivial infrared-safe prediction. For our problem in the classical soft limit, as we will see, only the 1-point event shape will be non-trivial because of factorization properties. Let's first discuss how to achieve the classical limit. We are interested in the limit where those gravitons contribute coherently to a (suitably localized) wave travelling in a specific direction. To make this precise, one can study the interaction hamiltonian for our theory of scalars minimally coupled with gravity in the linearized gravity approximation and in asymptotic limit |t| → +∞ [135] where H 0 is the free hamiltonian, is the number operator for the ith scalar field and (t, x) are the coordinates of the asymptotic trajectory of the scalar particles. This result is is valid only at first order in the perturbation h µν , but it tells something important about the asymptotic dynamics: solving the asymptotic Schrödinger equation for the potential V asy (t) gives an evolution operator which generates an infinite number of (soft) gravitons, or more precisely a coherent state of gravitons 25 for a suitable function α p (p) related to the structure of the asymptotic stress momentum tensor of the scalar particle T asy µν . In this specific context the state (7.24) is referred as Faddeev-Kulish state [136]. This gives an insight about how the quantum scattering theory generates a classical gravitational wave, composed of infinitely many gravitons, at least in the soft kinematic regime. More details about the connection of coherent states with the classical limit will be provided in [126].
In our case we can take the sum over all contributions M = 1, .., +∞, which implies to sum over M = 0, .., +∞ in our expression for the gravitational wave event shape in (7.17): . This is consistent with the graviton production rate in Weinberg [64] under our hypotheses. Once normalized by the total cross section we have There is also an additional term which is the analogue of the Coulomb phase contribution in QED, which plays an important role at higher orders in the soft expansion but not at leading order. and if we try to compute the gravitational energy-energy correlator at the classical level Ẽ GR (n 1 )Ẽ GR (n 2 )) class = Ẽ GR (n 1 ) class Ẽ GR (n 2 ) class , (7.27) which is a consequence of the classical (soft) factorization. This can be shown explicitly using the coherent state in (7.24).

Infrared safety of gravitational energy event shapes and soft BMS supertranslation
It is possible to do an interesting consistency check about infrared safety of the gravitational energy event shapes by looking at how soft supertranslation symmetry affects such matrix elements. Earlier we have used an infrared cutoff and we have observed the cancellation of infrared divergences as in the Bloch-Nordsieck mechanism, but it would be nice to understand if there is a theoretical reason for these event shapes to be infrared finite.
It is well known that the action of the soft supertranslation mode corresponds to an insertion of a soft graviton at the level of S-matrix elements [110,112]. In particular the generator of such symmetry reads where f (ζ,ζ) is an arbitrary function on the celestial sphere. Its action on graviton creation operators is given by 29) and the one with creation operators is fixed because T soft (f ) is hermitian. The latter implies that We now consider a gravitational event shape of the form  (7.32) and since the addition of a soft graviton gives the usual Weinberg soft factor with a pole in E p , the action (7.32) annihilates every S-matrix element (once we insert the completeness relation) [112].
It would be interesting to understand in more details if the invariance under the soft modes of large gauge symmetries, BMS supertranslations and in general asymptotic symmetries of massless particles always guarantees the IR-finiteness (in the soft sense 26 ) of a matrix element in a QFT. This could be helpful to prove in general whether for other onshell observables, for example like the ones defined in [57], infrared divergences are always going to cancel.

Discussion and Conclusions
Einstein equations solved near null infinity and BMS balance laws are of particular relevance for studying the property of the radiation at infinity, like its energy and its angular momentum which are captured by detectors localized at spatial infinity in a particular direction n. Once we integrate Einstein equations over the retarded time, say in Bondi gauge, new non-local operators appear which are commonly called light-ray operators. For our definition of light-ray operators, we integrate over the entire retarded time interval u ∈ R, which allows us to extract a suitable notion of density of energy and angular momentum flux produced during a scattering event. This is the natural framework where contributions coming from different massless particles combine together in a single expression.
Besides the direct physical applications, the algebra of such non-local light-ray operators at infinity is of intrinsic interest since it can be considered some sort of "detector algebra". In the conformal case the algebra was clarified in a beautiful work by Cordova and Shao [73] where they took advantage of the Hofman-Maldacena map [9] to move the light-sheet at infinity to a finite distance and they analized the structure of the algebra based on pure symmetry arguments. Therefore, it is of particular relevance to understand such algebra in the case of physical theories (possibly non-conformal) in 4 dimensions, which is the main purpose of this work.
Here we interpret the insertion of light-ray operators in a correlation function as the insertion of a physical detector which is capturing properties of the distribution of outgoing particles at null infinity and we work in perturbation theory. In particular, we focus on the case of massless particles with integer spin and we use the saddle point estimate for the leading contribution to the stress tensor components. This provides an explicit representation of the hard mode contributions to light-ray operators, which is of particular relevance for physical detectors which have any kind of lower energy resolution. The structure we find is consistent with the quantum fields localizing on point-particle trajectory, since the stress tensor expression reduces to the pointlike one, possibly with the spin contribution [99]. The universality of these expressions seems to point to some universality of light-ray operators in 4 dimensions beyond the conformal case. While a more rigorous analysis would be desired for the [N z (n 1 ), Nz(n 2 )] commutator, our results are consistent with the algebra of Cordova and Shao [73]. For the spin 1 case coupled with matter, an interaction term coming from the mixing of the so-called radiative and coulombic modes [101,102] seems to break the algebra. 26 Here we are not considering collinear divergences, which play a huge role in QCD.
There could be extensions of the algebra in the soft sector for the angular momentum flux, in particular in the non-abelian case where we have observed a non-trivial boundary contribution whose insertion is related to the soft gluon theorem. This would required a detailed analysis of the Hilbert space factorization into soft and hard modes, which is beyond the scope of the paper [93]. More generally it would be interesting to establish a basis of non-local operators at null infinity, making contact also with integrals of current operators along the retarded time and to Wilson line operators. It would be also nice to have a full analysis of the algebra in the Bondi gauge, for example by using Newman-Penrose techniques [137,138].
Borrowing some ideas from the QCD event shape detector operators, we have defined a suitable notion of energy event shape in the gravitational case using the Bondi news squared operator. Working with an EFT of gravity in the classical limit, we study the gravitational radiation produced in a particular directionn in the sky from the scattering of compact objects in the leading soft approximation. In this simple limit, the Glauber-Sudarshan representation [123,125] of the QFT radiation density matrix reduces to a single coherent state and we can sum up all graviton contributions to get the classical result which matches the literature on the zero frequency limit (ZFL) [61,70]. Working in the approximation where we neglect the single graviton energies compared to the energy of the wave, we observe that there is a cancellation of infrared divergences for the gravitational energy event shape which gives an infrared finite result. It would be interesting to extend our analysis to include classical corrections to our result due to additional finite energy contributions and outside the zero frequency limit, as done in [70,139] for the soft gravitational radiation. In the classical limit, the correlation function of gravitational energy event shapes factorize because the graviton emissions are uncorrelated, which is in agreement with the single coherent state structure. In particular we have defined also the 2-point gravitational energy correlator, which is particularly interesting especially outside the soft classical limit where non-trivial contributions should come into play. Since it should be a non-trivial infrared safe prediction of quantum gravity, it would be very fascinating to explore this further. Interestingly, recently it was pointed out that a signature of the quantum state of gravitons might be observable especially if it deviates from the coherent state prediction [140][141][142].

B Detector light-cone definition ofẼ(n) andK(n) light-ray operators
To connect our discussion to the usual QCD event shape literature and in order to indroduce a more covariant formalism, we identify a detector by its set of coordinates x µ = (t, rn) = (x − n µ + x +n µ ) where We can complete the basis with two spacelike vectors 27 m A which satisfy so that a generic point in spacetime is uniquely identified as x µ = x − n µ + x +n µ + x A m µ A . The completeness relation reads η µν = n µnν +n µ n ν (n ·n) which implies that x · y = 2(x − y + + x + y − ) − x A y A and that the flat metric can be written as It is worth pointing out that since we are interested in the algebra of different detectors it makes more sense to use a system of coordinates adapted at null infinity, as we did in the main text. Moreover, the latter extends to the more general case of asymptotically flat geometries and that is crucial for physical applications. We can write the ANEC operator at infinity in a covariant form as where we have defined a covariant light-cone projection T −− =n µnν (n·n) 2 T µν . This expression is Lorentz covariant under rescalings of the type n → ρ nn → ρn (B.5) and it has weight −3. It is important to stress that for external non-localized states, like momentum eigenstates, we need to take the large x − limit first and then perform the integral as done in (B.4). The opposite is usually done in CFT for localized states, for which all notions are equivalent [9,11]. In a similar way, we can write the momentum flow operator in a Lorentz covariant form as The next operator which we want to define isK(n) which has weight −2. The formulation of the operatorÑ A (n), where the index A has to be interpreted as a projection over the transverse direction m A (Ñ A (n) = (m A ) µÑ µ(n) (m A ·m A ) ), is more difficult. This is related to the fact that one needs to rescale the index A with x − to have a well-defined limit at infinity (see the related section 6.4).

C Light-ray operators for a light-sheet in the bulk
In the general, one would like to generalize the light-ray operator definition to allow them to be localized on general light-sheets in the bulk of the spacetime. If the theory is scale invariant (i.e. for all CFT for example) then the two definitions will coincide (see the explicit Hofman-Maldacena conformal transformation [9]) but in general they will be different. To simplify our calculations, we can definẽ so that light-sheets at fixedx − = const will be parametrized by the coordinatex + . It is worth noticing that asx − → +∞ we will havex ± → x ± , so that there is no ambiguity with the previous notation in appendix B. For the light-sheet defined byx − = c with c ∈ R, we define the light-ray family as The light-ray operators are then defined as It is quite remarkable that interactions drop out completely. This generically happens for scalar theories without derivative couplings and it implies some universality of the scalar light-ray algebra.

C.2 Spin 1
We have , where in the last line have imposed the light-cone gauge We see that while T −− (x) is interaction independent in this gauge, T A− (x) is not. Nevertheless, we are interested in light-ray operators where we integrate over the light-cone timẽ x + and therefore we can take advantage of the equations of motion provided that we take into account possible boundary terms atx + → ±∞. To fix completely the gauge for the boundary contributions, we can enforce the radiation gauge where together with the light-cone condition A − = 0 we have The light-ray operators are then defined as (C.10) If we compare to the self-interacting scalar case, we see that in YM theory interactions do not affect the definition of the operators E gluon (x ⊥ ) and K gluon (x ⊥ ) but they appear in N A,gluon (x ⊥ ) through the boundary contribution.

D Shear-inclusive ANEC and Raychaudhuri equation
The connection of the ANEC operator and its expectation value with causality properties of the underlying quantum field theory is well-known [144]. A congruence of null (complete, achronal) geodesics obeys Raychaudhuri equation where λ is the affine parameter, Θ = 1 A dA dλ 28 is the expansion, σ µν is the shear and k µ is a null vector of unit affine length which in our case we will take the light-sheet generator k µ = δ µ − . The right hand side of (D.1) immediately reminds us the enlarged definition of the matter ANEC with the gravitational contribution once we integrate along the null direction of the generator, albeit here (D.1) is purely classical in its nature. If we consider a generic congruence of null geodesics with starts at I − and ends on I + , one can ask how the expansion parameter changes along the congruence. Classically, the ANEC condition T matter µν k µ k ν ≥ 0 implies the classical focussing theorem, that is light-rays never "antifocus" as long as matter has positive energy. Quantum mechanically, quantum fields can developed negative energy so that we generically need a quantum version of the ANEC [145] and of the focussing theorem [146]. In general, there are problems in including quantized gravity perturbations in this story. 29 Nevertheless, assuming we have matter minimally coupled with gravity, for a stationary null surface one can quantize h µν on the null surface in light-front quantization [88] and at the first order in the perturbation one can define a shear-inclusive ANEC operator of the form E shear-inclusive = dλ 1 8πG σ µν σ µν + T matter µν k µ k ν , (D.2) whose expectation value is positive definite E shear-inclusive ≥ 0 as it can be proved from the linearized Raychaudhuri equation [147]. The connection to our previous definition of the shear-inclusive ANEC in the Bondi gauge consists in considering the asymptotically flat region for the geodesic congruence. Indeed in such region one can define a boundary shear as the limit of the rescaled shear of the congruence on future null infinity so that lim r→+∞ r 2 1 8πG σ µν σ µν + T matter µν k µ k ν ∼ 1 32πG N ζζ N ζζ + lim as proved in [148,149] in order to study asymptotic entropy bounds. 28 A is the cross sectional area here. 29 We thank Aron Wall for useful discussions on this point.