Skip to content

Commit

Permalink
Added includes and fixed missing M_PIf
Browse files Browse the repository at this point in the history
  • Loading branch information
EdoPro98 committed Oct 23, 2023
1 parent 6f2ac9a commit 67f615f
Show file tree
Hide file tree
Showing 5 changed files with 30 additions and 27 deletions.
1 change: 1 addition & 0 deletions include/SiPMAnalogSignal.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,7 @@
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <sstream>
#include <vector>

namespace sipm {
Expand Down
1 change: 1 addition & 0 deletions include/SiPMDebugInfo.h
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <sstream>

namespace sipm {
struct SiPMDebugInfo {
Expand Down
1 change: 1 addition & 0 deletions include/SiPMHit.h
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <sstream>

namespace sipm {
class SiPMHit {
Expand Down
1 change: 1 addition & 0 deletions include/SiPMProperties.h
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
#include <iomanip>
#include <iostream>
#include <map>
#include <sstream>
#include <string>
#include <vector>

Expand Down
53 changes: 26 additions & 27 deletions src/SiPMRandom.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,6 @@
#include "SiPMTypes.h"
#include <cstdint>


// Random seed
#include <x86intrin.h>
static uint64_t rdtsc() {
Expand Down Expand Up @@ -301,30 +300,30 @@ uint32_t SiPMRandom::randInteger(const uint32_t max) noexcept {
}

/**
* @param n Number of values to generate
*/
std::vector<double> SiPMRandom::Rand(const uint32_t n) {
std::vector<double> f64(n);
uint64_t* u64; // Free at end of function

if (n > 8) {
// simd8 allocates data
u64 = m_rng.simd8(n);
} else {
// We do allocation
u64 = (uint64_t*)aligned_alloc(64, n * sizeof(uint64_t));
for (uint32_t i = 0; i < n; ++i) {
u64[i] = m_rng();
}
}
* @param n Number of values to generate
*/
std::vector<double> SiPMRandom::Rand(const uint32_t n) {
std::vector<double> f64(n);
uint64_t* u64; // Free at end of function

if (n > 8) {
// simd8 allocates data
u64 = m_rng.simd8(n);
} else {
// We do allocation
u64 = (uint64_t*)aligned_alloc(64, n * sizeof(uint64_t));
for (uint32_t i = 0; i < n; ++i) {
f64[i] = static_cast<double>(u64[i]) * 0x1p-64;
u64[i] = m_rng();
}
}

free(u64);
for (uint32_t i = 0; i < n; ++i) {
f64[i] = static_cast<double>(u64[i]) * 0x1p-64;
}

return f64;
free(u64);

return f64;
}

/**
Expand All @@ -338,7 +337,7 @@ std::vector<float> SiPMRandom::RandF(const uint32_t n) {
// simd8 allocates data
// need n/2 uint64_t to generate n uint32_t
// odd case requires one more generation for tail
u32 = (uint32_t*)m_rng.simd8((n+1) / 2);
u32 = (uint32_t*)m_rng.simd8((n + 1) / 2);
} else {
// We do allocation
u32 = (uint32_t*)aligned_alloc(64, n * sizeof(uint32_t));
Expand All @@ -353,7 +352,7 @@ std::vector<float> SiPMRandom::RandF(const uint32_t n) {
// See SiPMRandom.h for details
f32[i] = static_cast<float>(u32[i]) * 0x1p-32;
}
f32[n-1] = Rand<float>();
f32[n - 1] = Rand<float>();

free(u32);

Expand Down Expand Up @@ -393,7 +392,7 @@ std::vector<double> SiPMRandom::randGaussian(const double mu, const double sigma
std::vector<float> SiPMRandom::randGaussianF(const float mu, const float sigma, const uint32_t n) {
std::vector<float> out(n);
const std::vector<float> u = RandF(n);
constexpr float TWO_PI = 2 * M_PIf;
constexpr float TWO_PI = 2 * M_PI;

for (uint32_t i = 0; i < n - 1; i += 2) {
const float R = sqrtf(-2.0 * logf(u[i]));
Expand All @@ -417,18 +416,18 @@ std::vector<uint32_t> SiPMRandom::randInteger(const uint32_t max, const uint32_t
std::vector<uint32_t> out(n);
uint32_t* u32;

if(n < 16){
if (n < 16) {
u32 = (uint32_t*)aligned_alloc(64, n * sizeof(uint32_t));
for(uint32_t i = 0; i < n; ++i){
for (uint32_t i = 0; i < n; ++i) {
u32[i] = m_rng() >> 32;
}
} else {
u32 = (uint32_t*)m_rng.simd8((n+1)/2);
u32 = (uint32_t*)m_rng.simd8((n + 1) / 2);
}

// Sort of fixed point arithmetic
// Avoids division and float numbers
for (uint32_t i = 0; i < n ; ++i) {
for (uint32_t i = 0; i < n; ++i) {
const uint64_t m = uint64_t(u32[i]) * uint64_t(max);
out[i] = m >> 32;
}
Expand Down

0 comments on commit 67f615f

Please sign in to comment.