Add units to StochasticRaytracing

src/ra/acoustic/StochasticRaytracing.cpp

@@ -7,17 +7,17 @@ namespace ra {
 
 StochasticRaytracing::StochasticRaytracing(Room room) noexcept : _room{std::move(room)} {}
 
-auto StochasticRaytracing::operator()(Simulation const& simulation) const -> Result
+auto StochasticRaytracing::operator()(Simulation const& sim) const -> Result
 {
     auto rng = std::mt19937{std::random_device{}()};
 
-    auto const rays         = randomRaysOnSphere(simulation.rays, rng);
-    auto const numTimeSteps = static_cast<std::size_t>(simulation.duration / simulation.timeStep);
+    auto const rays         = randomRaysOnSphere(sim.rays, rng);
+    auto const numTimeSteps = static_cast<std::size_t>((sim.duration / sim.timeStep).numerical_value_in(one));
 
-    auto histogram = std::vector(simulation.frequencies.size(), std::vector<double>(numTimeSteps));
+    auto histogram = std::vector(sim.frequencies.size(), std::vector<double>(numTimeSteps));
     for (auto frequency{0UL}; frequency < histogram.size(); ++frequency) {
         for (auto const& ray : rays) {
-            tarceRay(simulation, ray, histogram[frequency], frequency, rng);
+            tarceRay(sim, ray, histogram[frequency], frequency, rng);
         }
     }
 
@@ -45,7 +45,7 @@ auto StochasticRaytracing::tarceRay(
 
     // Initialize ray travel time. Ray tracing is terminated when the
     // travel time exceeds the impulse response length.
-    auto rayTime = std::chrono::duration<double>{0.0};
+    auto rayTime = 0.0 * si::second;
 
     // Initialize the ray energy to a normalized value of 1. Energy
     // decreases when the ray hits a surface.
@@ -67,7 +67,7 @@ auto StochasticRaytracing::tarceRay(
         rayPos = impactPosition;
 
         // Update cumulative ray travel time
-        rayTime = rayTime + std::chrono::duration<double>{distance / speedOfSound};
+        rayTime = rayTime + (distance / speedOfSound) * si::second;
 
         // Apply surface reflection to ray's energy
         // This is the amount of energy that is not lost through
@@ -84,7 +84,7 @@ auto StochasticRaytracing::tarceRay(
 
         // Determine the ray's time of arrival at receiver.
         auto const recvDistance  = std::sqrt(sum(recvDir * recvDir));
-        auto const timeOfArrival = rayTime + std::chrono::duration<double>{recvDistance / speedOfSound};
+        auto const timeOfArrival = rayTime + (recvDistance / speedOfSound) * si::second;
 
         if (timeOfArrival > sim.duration) {
             return;
@@ -93,15 +93,16 @@ auto StochasticRaytracing::tarceRay(
         // Determine amount of diffuse energy that reaches the receiver. See (5.20) in [2].
 
         // Compute received energy
-        auto const radius       = sim.radius;
+        auto const radius       = sim.radius.numerical_value_in(si::metre);
         auto const recvDir2     = pow(recvDir, 2.0);
         auto const impactNormal = getWallNormal(surfaceIdx);
         auto const cosTheta     = sum(recvDir * impactNormal) / (sqrt(sum(recvDir2)));
         auto const cosAlpha     = sqrt(sum(recvDir2) - std::pow(radius, 2.0)) / sum(recvDir2);
         auto const energy       = (1 - cosAlpha) * 2 * cosTheta * recvEnergy;
 
         // Update energy histogram
-        auto const timeIdx = static_cast<size_t>(std::max(0L, std::lround(timeOfArrival / sim.timeStep) - 1));
+        auto const idx     = std::lround((timeOfArrival / sim.timeStep).numerical_value_in(one));
+        auto const timeIdx = static_cast<size_t>(std::max(0L, idx - 1));
         histogram[timeIdx] = histogram[timeIdx] + energy;
 
         // Compute a new direction for the ray.

src/ra/acoustic/StochasticRaytracing.hpp

@@ -28,11 +28,11 @@ struct StochasticRaytracing
 
     struct Simulation
     {
-        std::vector<double> frequencies;
+        std::vector<quantity<isq::frequency[si::hertz]>> frequencies;
+        quantity<isq::duration[si::second]> duration;
+        quantity<isq::duration[si::second]> timeStep;
+        quantity<isq::radius[si::metre]> radius;
         std::size_t rays;
-        std::chrono::duration<double> duration;
-        std::chrono::duration<double> timeStep;
-        double radius;
     };
 
     using Result = std::vector<std::vector<double>>;

tool/RaumAkustik/tabs/stochastic_raytracing_tab.cpp

@@ -62,12 +62,27 @@ auto StochasticRaytracingEditor::resized() -> void
 
 auto StochasticRaytracingEditor::run() -> void
 {
+    using si::unit_symbols::Hz;
+
+    auto const frequencies = std::vector<quantity<isq::frequency[si::hertz]>>{
+        31.25 * Hz,
+        62.5 * Hz,
+        125.0 * Hz,
+        250.0 * Hz,
+        500.0 * Hz,
+        1000.0 * Hz,
+        2000.0 * Hz,
+        4000.0 * Hz,
+        8000.0 * Hz,
+        16'000.0 * Hz,
+    };
+
     auto const simulation = StochasticRaytracing::Simulation{
-        .frequencies = std::vector{31.25, 62.5, 125.0, 250.0, 500.0, 1000.0, 2000.0, 4000.0, 8000.0, 16'000.0},
+        .frequencies = frequencies,
+        .duration    = static_cast<double>(_duration.getValue()) * si::second,
+        .timeStep    = 0.001 * si::second,
+        .radius      = 0.0875 * si::metre,
         .rays        = static_cast<size_t>(static_cast<double>(_rays.getValue())),
-        .duration    = std::chrono::duration<double>{static_cast<double>(_duration.getValue())},
-        .timeStep    = std::chrono::duration<double>{0.001},
-        .radius      = 0.0875,
     };
 
     auto const roomLayout      = _roomEditor.getRoomLayout();
@@ -111,8 +126,8 @@ auto StochasticRaytracingEditor::run() -> void
     }
 
     for (auto i{0U}; i < _result->size(); ++i) {
-        _plots[static_cast<int>(i)]
-            ->plot(juce::String(simulation.frequencies[i]) + "Hz", _result->at(i), simulation.duration, _maxGain);
+        auto const f = simulation.frequencies[i].numerical_value_in(si::hertz);
+        _plots[static_cast<int>(i)]->plot(juce::String(f) + "Hz", _result->at(i), simulation.duration, _maxGain);
     }
 
     repaint();