Skip to main content

laddu_core/variables/
polarization.rs

1use std::fmt::Display;
2
3use serde::{Deserialize, Serialize};
4
5use super::{AuxSelection, Variable};
6use crate::{
7    data::{DatasetMetadata, EventLike},
8    reaction::Reaction,
9    vectors::Vec3,
10    LadduResult,
11};
12
13/// A struct defining the polarization angle for a beam relative to the production plane.
14#[derive(Clone, Debug, Serialize, Deserialize)]
15pub struct PolAngle {
16    reaction: Reaction,
17    angle_aux: AuxSelection,
18}
19
20impl Display for PolAngle {
21    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
22        write!(
23            f,
24            "PolAngle(reaction={:?}, angle_aux={})",
25            self.reaction.topology(),
26            self.angle_aux.name(),
27        )
28    }
29}
30
31impl PolAngle {
32    /// Constructs the polarization angle given a [`Reaction`] describing the production plane and
33    /// the auxiliary column storing the precomputed angle.
34    pub fn new<A>(reaction: Reaction, angle_aux: A) -> Self
35    where
36        A: Into<String>,
37    {
38        Self {
39            reaction,
40            angle_aux: AuxSelection::new(angle_aux.into()),
41        }
42    }
43}
44
45#[typetag::serde]
46impl Variable for PolAngle {
47    fn bind(&mut self, metadata: &DatasetMetadata) -> LadduResult<()> {
48        let _ = metadata;
49        self.angle_aux.bind(metadata)?;
50        Ok(())
51    }
52
53    fn value(&self, event: &dyn EventLike) -> f64 {
54        let resolved = self
55            .reaction
56            .resolve_two_to_two(event)
57            .unwrap_or_else(|err| panic!("failed to evaluate polarization angle: {err}"));
58        let beam = resolved.p1();
59        let recoil = resolved.p4();
60        let pol_angle = event.aux_at(self.angle_aux.index());
61        let polarization = Vec3::new(pol_angle.cos(), pol_angle.sin(), 0.0);
62        let y = beam.vec3().cross(&-recoil.vec3()).unit();
63        let numerator = y.dot(&polarization);
64        let denominator = beam.vec3().unit().dot(&polarization.cross(&y));
65        f64::atan2(numerator, denominator)
66    }
67}
68
69/// A struct defining the polarization magnitude for a beam relative to the production plane.
70#[derive(Clone, Debug, Serialize, Deserialize)]
71pub struct PolMagnitude {
72    magnitude_aux: AuxSelection,
73}
74
75impl Display for PolMagnitude {
76    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
77        write!(
78            f,
79            "PolMagnitude(magnitude_aux={})",
80            self.magnitude_aux.name(),
81        )
82    }
83}
84
85impl PolMagnitude {
86    /// Constructs the polarization magnitude given the named auxiliary column containing the
87    /// magnitude value.
88    pub fn new<S: Into<String>>(magnitude_aux: S) -> Self {
89        Self {
90            magnitude_aux: AuxSelection::new(magnitude_aux.into()),
91        }
92    }
93}
94
95#[typetag::serde]
96impl Variable for PolMagnitude {
97    fn bind(&mut self, metadata: &DatasetMetadata) -> LadduResult<()> {
98        self.magnitude_aux.bind(metadata)
99    }
100
101    fn value(&self, event: &dyn EventLike) -> f64 {
102        event.aux_at(self.magnitude_aux.index())
103    }
104}
105
106/// A struct for obtaining both the polarization angle and magnitude at the same time.
107#[derive(Clone, Debug, Serialize, Deserialize)]
108pub struct Polarization {
109    /// See [`PolMagnitude`].
110    pub pol_magnitude: PolMagnitude,
111    /// See [`PolAngle`].
112    pub pol_angle: PolAngle,
113}
114
115impl Display for Polarization {
116    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
117        write!(
118            f,
119            "Polarization(reaction={:?}, magnitude_aux={}, angle_aux={})",
120            self.pol_angle.reaction.topology(),
121            self.pol_magnitude.magnitude_aux.name(),
122            self.pol_angle.angle_aux.name(),
123        )
124    }
125}
126
127impl Polarization {
128    /// Constructs the polarization angle and magnitude given a [`Reaction`] and distinct
129    /// auxiliary columns for magnitude and angle.
130    ///
131    /// # Panics
132    ///
133    /// Panics if `magnitude_aux` and `angle_aux` refer to the same auxiliary column name.
134    pub fn new<M, A>(reaction: Reaction, magnitude_aux: M, angle_aux: A) -> Self
135    where
136        M: Into<String>,
137        A: Into<String>,
138    {
139        let magnitude_aux = magnitude_aux.into();
140        let angle_aux = angle_aux.into();
141        assert!(
142            magnitude_aux != angle_aux,
143            "Polarization magnitude and angle must reference distinct auxiliary columns"
144        );
145        Self {
146            pol_magnitude: PolMagnitude::new(magnitude_aux),
147            pol_angle: PolAngle::new(reaction, angle_aux),
148        }
149    }
150}