Commit 07fc1a1b authored by Roel Aaij's avatar Roel Aaij
Browse files

Add some preliminary tests.

parent e7957243
...@@ -71,9 +71,8 @@ std::tuple<storage_t, storage_t> generate_avx2(const long time_start, const long ...@@ -71,9 +71,8 @@ std::tuple<storage_t, storage_t> generate_avx2(const long time_start, const long
const float tau_l0 = gens.tau_l0(); const float tau_l0 = gens.tau_l0();
size_t storage_size = n_expect + long_v::size() - n_expect % long_v::size(); size_t storage_size = n_expect + long_v::size() - n_expect % long_v::size();
storage_t times; times.resize(storage_size); storage_t times; times.resize(n_expect);
storage_t values; values.resize(storage_size); storage_t values; values.resize(n_expect + 2 * long_v::size());
times[0] = 0l;
auto int_v_to_long_v = [] (const int_v& in) -> pair<long_v, long_v> auto int_v_to_long_v = [] (const int_v& in) -> pair<long_v, long_v>
{ {
...@@ -145,13 +144,13 @@ std::tuple<storage_t, storage_t> generate_avx2(const long time_start, const long ...@@ -145,13 +144,13 @@ std::tuple<storage_t, storage_t> generate_avx2(const long time_start, const long
auto val0 = simd_cast<int_v>(z0) | (pmt1 << 8) | ((100 * (dom + 1) + mod + 1) << 13); auto val0 = simd_cast<int_v>(z0) | (pmt1 << 8) | ((100 * (dom + 1) + mod + 1) << 13);
auto [val0_first, val0_second] = int_v_to_long_v(val0); auto [val0_first, val0_second] = int_v_to_long_v(val0);
val0_first.store(&values[idx - 2 * long_v::size()]); val0_first.store(&values[vidx]);
val0_second.store(&values[idx - long_v::size()]); val0_second.store(&values[vidx + long_v::size()]);
auto val1 = simd_cast<int_v>(z1) | (pmt2 << 8) | ((100 * (dom + 1) + mod + 1) << 13); auto val1 = simd_cast<int_v>(z1) | (pmt2 << 8) | ((100 * (dom + 1) + mod + 1) << 13);
auto [val1_first, val1_second] = int_v_to_long_v(val1); auto [val1_first, val1_second] = int_v_to_long_v(val1);
val1_first.store(&values[idx]); val1_first.store(&values[vidx + 2 * long_v::size()]);
val1_second.store(&values[idx + long_v::size()]); val1_second.store(&values[vidx + 3 *long_v::size()]);
} }
} }
} }
......
...@@ -50,8 +50,7 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo ...@@ -50,8 +50,7 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo
const float tau_l0 = gens.tau_l0(); const float tau_l0 = gens.tau_l0();
storage_t times; times.resize(n_expect); storage_t times; times.resize(n_expect);
storage_t values; values.resize(n_expect); storage_t values; values.resize(n_expect + 1);
times[0] = 0l;
size_t idx = 0; size_t idx = 0;
...@@ -73,9 +72,7 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo ...@@ -73,9 +72,7 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo
fill_coincidences(times, idx, time_start, time_end, gens); fill_coincidences(times, idx, time_start, time_end, gens);
// fill values // fill values
const size_t value_end = idx + (2 - (idx % 2)); for (size_t vidx = mod_start; vidx < idx; vidx += 2) {
for (size_t vidx = mod_start; vidx < value_end; vidx += 2) {
long pmt1 = flat_pmt(mt); long pmt1 = flat_pmt(mt);
long pmt2 = flat_pmt(mt); long pmt2 = flat_pmt(mt);
...@@ -90,10 +87,10 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo ...@@ -90,10 +87,10 @@ std::tuple<storage_t, storage_t> generate_scalar(const long time_start, const lo
z1 = fact * tot_sigma * z1 + tot_mean; z1 = fact * tot_sigma * z1 + tot_mean;
auto val0 = static_cast<int>(z0) | (pmt1 << 8) | ((100 * (dom + 1) + mod + 1) << 13); auto val0 = static_cast<int>(z0) | (pmt1 << 8) | ((100 * (dom + 1) + mod + 1) << 13);
values[idx - 1] = val0; values[vidx] = val0;
auto val1 = static_cast<int>(z1) | (pmt2 << 8) | ((100 * (dom + 1) + mod + 1) << 13); auto val1 = static_cast<int>(z1) | (pmt2 << 8) | ((100 * (dom + 1) + mod + 1) << 13);
values[idx] = val1; values[vidx + 1] = val1;
} }
} }
} }
......
...@@ -23,6 +23,10 @@ endif() ...@@ -23,6 +23,10 @@ endif()
include_directories(include) include_directories(include)
if (USE_AVX2)
add_definitions("-DUSE_AVX2")
endif()
add_library(test_functions STATIC add_library(test_functions STATIC
lib/coincidences.cpp lib/coincidences.cpp
lib/generate.cpp) lib/generate.cpp)
...@@ -34,6 +38,16 @@ add_executable(test_k40gen ...@@ -34,6 +38,16 @@ add_executable(test_k40gen
target_link_libraries(test_k40gen PRIVATE generate test_functions) target_link_libraries(test_k40gen PRIVATE generate test_functions)
add_test(TestCoincidences test_coincidences) add_test(TestK40 test_k40gen)
find_package(ROOT COMPONENTS Core Cling Hist Gpad)
add_executable(test_k40gen_root
src/test_rates_root.cpp)
target_include_directories(test_k40gen_root PRIVATE
${ROOT_INCLUDE_DIRS})
find_package(ROOT) target_link_libraries(test_k40gen_root PRIVATE
generate
test_functions
${ROOT_LIBRARIES})
...@@ -2,6 +2,15 @@ ...@@ -2,6 +2,15 @@
#include <array> #include <array>
std::pair<double, double> generate_times(std::array<float, 4> rates, bool use_avx2); template <std::size_t N>
struct get_n {
template <typename T>
auto operator()(T&& t) const ->
decltype(std::get<N>(std::forward<T>(t))) {
return std::get<N>(std::forward<T>(t));
}
};
std::pair<double, double> generate_l0(float l0_rate, bool use_avx2);
std::pair<double, double> coincidence_rate(std::array<float, 4> rates); std::pair<double, double> coincidence_rate(std::array<float, 4> rates);
...@@ -4,19 +4,13 @@ ...@@ -4,19 +4,13 @@
#include <range/v3/view/zip.hpp> #include <range/v3/view/zip.hpp>
#include <range/v3/algorithm/sort.hpp> #include <range/v3/algorithm/sort.hpp>
template <std::size_t N> #include <test_functions.h>
struct get_n {
template <typename T>
auto operator()(T&& t) const ->
decltype(std::get<N>(std::forward<T>(t))) {
return std::get<N>(std::forward<T>(t));
}
};
using namespace std; using namespace std;
pair<double, double> generate_times(array<float, 4> rates, bool use_avx2) { pair<double, double> generate_l0(float l0_rate, bool use_avx2) {
array<float, 4> rates = {l0_rate, 0., 0., 0.};
Generators gens{1052, 9523, rates}; Generators gens{1052, 9523, rates};
long dt = std::lround(1e8); long dt = std::lround(1e8);
...@@ -35,7 +29,7 @@ pair<double, double> generate_times(array<float, 4> rates, bool use_avx2) { ...@@ -35,7 +29,7 @@ pair<double, double> generate_times(array<float, 4> rates, bool use_avx2) {
double av = 0.; double av = 0.;
const size_t n_times = times.size(); const size_t n_times = times.size();
for (size_t i = 0; i < n_times - 1; ++i) { for (size_t i = 0; i < n_times - 1; ++i) {
av += static_cast<double>(times[i]) / n_times; av += static_cast<double>(times[i + 1] - times[i]) / n_times;
} }
return {av, times.size()}; return {av, times.size()};
......
...@@ -8,14 +8,19 @@ ...@@ -8,14 +8,19 @@
using namespace std; using namespace std;
TEST_CASE( "Rates make sense", "[rates]" ) { TEST_CASE( "Rates make sense", "[rates]" ) {
auto check_rates = [](array<float, 4> rates, bool use_avx2) { auto check_rates = [](float l0_rate, bool use_avx2) {
const auto [av_rate, av_n] = generate_times(rates, use_avx2); const auto [av_rate, av_n] = generate_l0(l0_rate, use_avx2);
auto coincidence_rate = std::accumulate(begin(rates), end(rates), 0.); REQUIRE(std::abs(l0_rate - av_rate) / l0_rate < 1e-4);
REQUIRE(std::abs(coincidence_rate - av_rate) / coincidence_rate < 1e-4);
cout << av_rate << " " << av_n << endl; cout << av_rate << " " << av_n << endl;
}; };
check_rates({7000.}, false); for (auto [rate, avx2] : {make_pair(7000.f, false),
check_rates({1000., 100.}, false); make_pair(1000.f, false)
check_rates({7000., 700., 70., 7.}, false); #ifdef USE_AVX2
, make_pair(7000.f, true),
make_pair(1000.f, true)
#endif
}) {
check_rates(rate, avx2);
}
} }
#include <numeric>
#include <iostream>
#include <array>
#include <generate.h>
#include <storage.h>
#include <test_functions.h>
#include <TApplication.h>
#include <TH1.h>
#include <TCanvas.h>
using namespace std;
int main(int argc, char* argv[]) {
TApplication app{"test_rate", &argc, argv};
map<bool, unique_ptr<TH1I>> histos;
array<float, 4> rates{7000., 0., 0., 0.};
for (bool use_avx2 : {false, true}) {
string name = string{"diffs_"} + (use_avx2 ? "avx2" : "scalar");
auto r = histos.emplace(use_avx2, make_unique<TH1I>(name.c_str(), name.c_str(), 100, 0, 1000000));
auto& time_histo = r.first->second;
Generators gens{1052, 9523, rates};
long dt = std::lround(1e8);
long time_start = 0, time_end = time_start + dt;
auto [times, values] = generate(time_start, time_end, gens, use_avx2);
const size_t n_times = times.size();
for (size_t i = 0; i < n_times - 1; ++i) {
if (((values[i + 1]) >> 13) == (values[i] >> 13)) {
time_histo->Fill(times[i + 1] - times[i]);
}
}
}
TCanvas canvas{"canvas", "canvas", 600, 800};
canvas.Divide(1, 2);
for (auto& [arch, histo] : histos) {
canvas.cd(arch + 1);
histo->Draw();
}
app.Run();
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment