autd3_core/acoustics/
mod.rs

1/// directivity module
2pub mod directivity;
3
4use std::f32::consts::PI;
5
6use crate::{
7    defined::T4010A1_AMPLITUDE,
8    geometry::{Complex, Point3, Transducer, UnitVector3},
9};
10
11use directivity::Directivity;
12
13/// Calculate the pressure at the target position.
14#[inline]
15#[must_use]
16pub fn propagate<D: Directivity>(
17    tr: &Transducer,
18    wavenumber: f32,
19    dir: &UnitVector3,
20    target_pos: &Point3,
21) -> Complex {
22    const P0: f32 = T4010A1_AMPLITUDE / (4. * PI);
23    let diff = target_pos - tr.position();
24    let dist = diff.norm();
25    Complex::from_polar(
26        P0 / dist * D::directivity_from_dir(dir, &diff),
27        wavenumber * dist,
28    )
29}
30
31#[cfg(test)]
32mod tests {
33    use super::*;
34
35    use rand::Rng;
36
37    use crate::{
38        defined::mm,
39        geometry::{Device, UnitQuaternion, Vector3},
40    };
41    use directivity::tests::TestDirectivity;
42
43    macro_rules! assert_complex_approx_eq {
44        ($a:expr, $b:expr) => {
45            approx::assert_abs_diff_eq!($a.re, $b.re, epsilon = 1e-3 / mm);
46            approx::assert_abs_diff_eq!($a.im, $b.im, epsilon = 1e-3 / mm);
47        };
48    }
49
50    #[rstest::fixture]
51    fn tr() -> Transducer {
52        let mut rng = rand::rng();
53        Transducer::new(
54            0,
55            0,
56            Point3::new(
57                rng.random_range(-100.0..100.0),
58                rng.random_range(-100.0..100.0),
59                rng.random_range(-100.0..100.0),
60            ),
61        )
62    }
63
64    #[rstest::fixture]
65    fn rot() -> UnitQuaternion {
66        let mut rng = rand::rng();
67        UnitQuaternion::from_axis_angle(
68            &Vector3::x_axis(),
69            rng.random_range::<f32, _>(-180.0..180.0).to_radians(),
70        ) * UnitQuaternion::from_axis_angle(
71            &Vector3::y_axis(),
72            rng.random_range::<f32, _>(-180.0..180.0).to_radians(),
73        ) * UnitQuaternion::from_axis_angle(
74            &Vector3::z_axis(),
75            rng.random_range::<f32, _>(-180.0..180.0).to_radians(),
76        )
77    }
78
79    #[rstest::fixture]
80    fn target() -> Point3 {
81        let mut rng = rand::rng();
82        Point3::new(
83            rng.random_range(-100.0..100.0),
84            rng.random_range(-100.0..100.0),
85            rng.random_range(-100.0..100.0),
86        )
87    }
88
89    #[rstest::fixture]
90    fn sound_speed() -> f32 {
91        let mut rng = rand::rng();
92        rng.random_range(300e3..400e3)
93    }
94
95    #[rstest::rstest]
96    #[test]
97    fn test_propagate(tr: Transducer, rot: UnitQuaternion, target: Point3, sound_speed: f32) {
98        let mut device = Device::new(0, rot, vec![tr.clone()]);
99        device.sound_speed = sound_speed;
100        let wavenumber = device.wavenumber();
101        assert_complex_approx_eq!(
102            {
103                let diff = target - tr.position();
104                let dist = diff.norm();
105                let r = T4010A1_AMPLITUDE / (4. * PI) / dist
106                    * TestDirectivity::directivity_from_dir(device.axial_direction(), &diff);
107                let phase = wavenumber * dist;
108                Complex::new(r * phase.cos(), r * phase.sin())
109            },
110            super::propagate::<TestDirectivity>(&tr, wavenumber, device.axial_direction(), &target)
111        );
112    }
113}