Skip to content

Commit

Permalink
feat: add ThreePiPolFrac amplitude (mirroring ThreePiAngles.cc in…
Browse files Browse the repository at this point in the history
… `halld_sim`)
  • Loading branch information
denehoffman committed Jul 31, 2024
1 parent f44585b commit 2e59d2f
Show file tree
Hide file tree
Showing 9 changed files with 355 additions and 7 deletions.
46 changes: 40 additions & 6 deletions Cargo.lock

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

1 change: 1 addition & 0 deletions Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ dyn-clone = "1.0.17"
tracing = "0.1.40"
ganesh = "0.4.0"
parking_lot = "0.12.3"
wigners = "0.3.0"

[profile.release]
lto = true
Expand Down
1 change: 1 addition & 0 deletions crates/rustitude-gluex/Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ pyo3 = { workspace = true }
rayon = { workspace = true }
rustitude-core = { workspace = true }
sphrs = { workspace = true }
wigners = { workspace = true }

[features]
default = []
Expand Down
2 changes: 1 addition & 1 deletion crates/rustitude-gluex/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@ There are a few amplitudes from [halld_sim](https://github.com/JeffersonLab/hall
| [PiPlusRegge.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/PiPlusRegge.cc) | | :x: |
| [Piecewise.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/Piecewise.cc) | `rustitude::amplitude::PiecewiseM` | :white_check_mark: |
| [SinglePS.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/SinglePS.cc) | `rustitude-gluex::harmonics::OnePS` | :white_check_mark: |
| [ThreePiAngles.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/ThreePiAngles.cc) | | :bangbang: |
| [ThreePiAngles.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/ThreePiAngles.cc) | `rustitude-gluex::polarization::ThreePiPolFrac` | :white_check_mark: |
| [ThreePiAnglesSchilling.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/ThreePiAnglesSchilling.cc) | `rustitude-gluex::sdmes::ThreePiSDME` | :white_check_mark: |
| [TwoLeptonAngles.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/TwoLeptonAngles.cc) | | :heavy_exclamation_mark: |
| [TwoLeptonAnglesGJ.cc](https://github.com/JeffersonLab/halld_sim/blob/master/src/libraries/AMPTOOLS_AMPS/TwoLeptonAnglesGJ.cc) | | :heavy_exclamation_mark: |
Expand Down
1 change: 1 addition & 0 deletions crates/rustitude-gluex/src/lib.rs
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
pub mod dalitz;
pub mod harmonics;
pub mod polarization;
pub mod resonances;
pub mod sdmes;
pub mod utils;
155 changes: 155 additions & 0 deletions crates/rustitude-gluex/src/polarization.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,155 @@
use crate::utils::{self, Decay, Frame, Sign};
use rayon::prelude::*;
use rustitude_core::prelude::*;
use sphrs::{ComplexSH, SHEval};

#[derive(Clone)]
pub struct ThreePiPolFrac<F> {
beam_pol: F,
j_resonance: u32,
p_resonance: F,
i_resonance: usize,
l_resonance: u32,
j_isobar: u32,
i_isobar: usize,
iz_daughters: [usize; 3],
decay_resonance: Decay,
decay_isobar: Decay,
data: Vec<Complex<F>>,
}

impl<F: Field> ThreePiPolFrac<F> {
#[allow(clippy::too_many_arguments)]
pub fn new(
beam_pol: Sign,
j_resonance: u32,
p_resonance: Sign,
i_resonance: usize,
l_resonance: u32,
j_isobar: u32,
i_isobar: usize,
iz_daughters: [usize; 3],
decay_resonance: Decay,
decay_isobar: Decay,
) -> Self {
match (decay_resonance, decay_isobar) {
(Decay::ThreeBodyDecay(_), Decay::TwoBodyDecay(_)) => Self {
beam_pol: match beam_pol {
Sign::Positive => F::ONE,
Sign::Negative => -F::ONE,
},
j_resonance,
p_resonance: match p_resonance {
Sign::Positive => F::ONE,
Sign::Negative => -F::ONE,
},
i_resonance,
l_resonance,
j_isobar,
i_isobar,
iz_daughters,
decay_resonance,
decay_isobar,
data: Vec::default(),
},
_ => unimplemented!(),
}
}
}

impl<F: Field> Node<F> for ThreePiPolFrac<F> {
fn calculate(&self, parameters: &[F], event: &Event<F>) -> Result<Complex<F>, RustitudeError> {
Ok(self.data[event.index] * (F::ONE + self.beam_pol * parameters[0]) / F::FOUR)
}

fn precalculate(&mut self, dataset: &Dataset<F>) -> Result<(), RustitudeError> {
self.data = dataset
.events
.par_iter()
.map(|event| {
let isobar_p4 = self.decay_isobar.resonance_p4(event);
let resonance_p4 = self.decay_resonance.resonance_p4(event);
let alpha = event.recoil_p4.phi();
let (_, _, _, p3_res_coords) =
self.decay_resonance.coordinates(Frame::Helicity, 2, event);
let p1_iso_p4 = self.decay_isobar.primary_p4(event).boost_along(&isobar_p4);
let (_, _, _, p1_iso_coords) =
Frame::Helicity.coordinates(self.decay_resonance, &p1_iso_p4, event);
let k = utils::breakup_momentum(
resonance_p4.m(),
isobar_p4.m(),
self.decay_resonance.tertiary_p4(event).m(),
);
let q = utils::breakup_momentum(
isobar_p4.m(),
self.decay_isobar.primary_p4(event).m(),
self.decay_isobar.secondary_p4(event).m(),
);
let neg_res_hel_prod =
Complex::new(F::fcos(F::TWO * alpha), F::fsin(F::TWO * alpha))
* self.beam_pol
* self.p_resonance
* (F::convert_u32((self.j_resonance % 2) * 2) / F::TWO);
let mut res = F::ZERO.c();
for m_res in -(self.l_resonance as isize)..(self.l_resonance as isize) {
let mut term = F::ZERO.c();
for m_iso in -(self.j_isobar as isize)..(self.j_isobar as isize) {
let ylm = ComplexSH::Spherical.eval(
self.j_isobar as i64,
m_iso as i64,
&p1_iso_coords,
);
let cg_neg = F::f(wigners::clebsch_gordan(
self.j_isobar,
self.l_resonance as i32,
m_iso as u32,
m_res as i32,
self.j_resonance,
-1,
));
let cg_pos = F::f(wigners::clebsch_gordan(
self.j_isobar,
self.l_resonance as i32,
m_iso as u32,
m_res as i32,
self.j_resonance,
1,
));
term += ylm * (neg_res_hel_prod * cg_neg + cg_pos);
}
let ylm = ComplexSH::Spherical.eval(
self.l_resonance as i64,
m_res as i64,
&p3_res_coords,
);
term *= ylm;
res += term;
}
res *= F::f(
wigners::clebsch_gordan(
1,
1,
self.iz_daughters[0] as u32,
self.iz_daughters[1] as i32,
self.i_isobar as u32,
(self.iz_daughters[0] + self.iz_daughters[1]) as i32,
) * wigners::clebsch_gordan(
self.i_isobar as u32,
1,
(self.iz_daughters[0] + self.iz_daughters[1]) as u32,
self.iz_daughters[2] as i32,
self.i_resonance as u32,
(self.iz_daughters[0] + self.iz_daughters[1] + self.iz_daughters[2]) as i32,
),
) * k.fpowi(self.l_resonance as i32)
* q.fpowi(self.j_isobar as i32);
res
})
.collect();
Ok(())
}

fn parameters(&self) -> Vec<String> {
std::vec!["polarization fraction".to_string()]
}
}
42 changes: 42 additions & 0 deletions py-rustitude/rustitude/gluex/polarization.pyi
Original file line number Diff line number Diff line change
@@ -0,0 +1,42 @@
from rustitude import Amplitude, Amplitude_64, Amplitude_32
from rustitude.gluex import Decay, Sign

def ThreePiPolFrac(
name: str,
beam_pol: Sign,
j_resonance: int,
p_resonance: int,
i_resonance: int,
l_resonance: int,
j_isobar: int,
i_isobar: int,
iz_daughters: tuple[int, int, int],
decay_resonance=Decay.ThreeBodyDecay([0, 1, 2]),
decay_isobar=Decay.TwoBodyDecay([0, 1]),
) -> Amplitude: ...
def ThreePiPolFrac_64(
name: str,
beam_pol: Sign,
j_resonance: int,
p_resonance: int,
i_resonance: int,
l_resonance: int,
j_isobar: int,
i_isobar: int,
iz_daughters: tuple[int, int, int],
decay_resonance=Decay.ThreeBodyDecay([0, 1, 2]),
decay_isobar=Decay.TwoBodyDecay([0, 1]),
) -> Amplitude_64: ...
def ThreePiPolFrac_32(
name: str,
beam_pol: Sign,
j_resonance: int,
p_resonance: int,
i_resonance: int,
l_resonance: int,
j_isobar: int,
i_isobar: int,
iz_daughters: tuple[int, int, int],
decay_resonance=Decay.ThreeBodyDecay([0, 1, 2]),
decay_isobar=Decay.TwoBodyDecay([0, 1]),
) -> Amplitude_32: ...
2 changes: 2 additions & 0 deletions py-rustitude/src/gluex/mod.rs
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
use pyo3::prelude::*;
mod dalitz;
mod harmonics;
mod polarization;
mod resonances;
mod sdmes;
use crate::add_submodule;
Expand All @@ -10,6 +11,7 @@ pub fn pyo3_module(m: &Bound<'_, PyModule>) -> PyResult<()> {
add_submodule(m, "rustitude.gluex.resonances", resonances::pyo3_module)?;
add_submodule(m, "rustitude.gluex.harmonics", harmonics::pyo3_module)?;
add_submodule(m, "rustitude.gluex.dalitz", dalitz::pyo3_module)?;
add_submodule(m, "rustitude.gluex.polarization", polarization::pyo3_module)?;
m.add_class::<rustitude_gluex::utils::Wave>()?;
m.add_class::<rustitude_gluex::utils::Frame>()?;
m.add_class::<rustitude_gluex::utils::Sign>()?;
Expand Down
Loading

0 comments on commit 2e59d2f

Please sign in to comment.