From 8bb10bbafa4402e9bb21369fa4cc32a1bcb824c1 Mon Sep 17 00:00:00 2001 From: Julien Date: Thu, 22 Feb 2024 16:51:40 +0100 Subject: [PATCH 1/4] add factor of $ into radial sampling --- src/Source.cpp | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/src/Source.cpp b/src/Source.cpp index 2ac28a100..9ff0b7f4d 100644 --- a/src/Source.cpp +++ b/src/Source.cpp @@ -362,7 +362,8 @@ void SourceUniformCylinder::setDescription() { // --------------------------------------------------------------------------- SourceSNRDistribution::SourceSNRDistribution() : - rEarth(8.5 * kpc),alpha(2.), beta(3.53), zg(0.3 * kpc) { + rEarth(8.5 * kpc), beta(3.53), zg(0.3 * kpc) { + setAlpha(2.); setFrMax(); setFzMax(0.3 * kpc); setRMax(20 * kpc); @@ -370,7 +371,8 @@ SourceSNRDistribution::SourceSNRDistribution() : } SourceSNRDistribution::SourceSNRDistribution(double rEarth, double alpha, double beta, double zg) : - rEarth(rEarth),alpha(alpha), beta(beta), zg(zg) { + rEarth(rEarth), beta(beta), zg(zg) { + setAlpha(alpha); setFrMax(); setFzMax(zg); setRMax(20 * kpc); @@ -414,7 +416,7 @@ double SourceSNRDistribution::fz(double z) const{ } double SourceSNRDistribution::getAlpha() const { - return alpha; + return alpha - 1; } double SourceSNRDistribution::getBeta() const { @@ -458,7 +460,7 @@ double SourceSNRDistribution::getZMax() const { } void SourceSNRDistribution::setAlpha(double a) { - alpha = a; + alpha = a + 1.; // add 1 for dV = r * dR * dphi * dz setRMax(rMax); setFrMax(); } @@ -504,7 +506,7 @@ void SourcePulsarDistribution::prepareParticle(ParticleState& particle) const { double Rtilde; while (true) { Rtilde = random.rand() * rMax; - double fTest = random.rand() * frMax; + double fTest = random.rand() * frMax * 1.1; double fR = fr(Rtilde); if (fTest <= fR) { break; @@ -530,10 +532,8 @@ void SourcePulsarDistribution::prepareParticle(ParticleState& particle) const { } double SourcePulsarDistribution::fr(double r) const { - double Atilde = (pow(beta, 4.) * exp(-beta)) / (12 * M_PI * pow(rEarth, 2.)); - double f = pow(r / rEarth, 2.) * exp(-beta * (r - rEarth) / rEarth); - double fr = Atilde * f; - return fr; + double f = r * pow(r / rEarth, 2.) * exp(-beta * (r - rEarth) / rEarth); + return f; } double SourcePulsarDistribution::fz(double z) const{ @@ -569,8 +569,8 @@ double SourcePulsarDistribution::blurTheta(double thetaTilde, double rTilde) con } void SourcePulsarDistribution::setFrMax(double R, double b) { - frMax = pow(b, 2.) / (3 * pow(R, 2.) * M_PI) * exp(-2.); - return; + double r = 3 * R / b; + frMax = fr(r); } void SourcePulsarDistribution::setFzMax(double zg) { From 096dafcfd9338c3a26815caa0a2fecaab0de7fa6 Mon Sep 17 00:00:00 2001 From: Julien Date: Thu, 14 Mar 2024 16:49:45 +0100 Subject: [PATCH 2/4] add bug-fix to changelog --- CHANGELOG.md | 1 + 1 file changed, 1 insertion(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index bf972ce31..ceabd8b5b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,7 @@ ### Bug fixes: * Fixed sign for exponential decay of magn. field strength with Galactic height in LogarithmicSpiralField + * Fixed r term in source distribution for SNR and Pulsar ### New features: From 0169c104cc730d336e61dc32f77dcdd51fe6f970 Mon Sep 17 00:00:00 2001 From: Julien Date: Thu, 14 Mar 2024 17:08:48 +0100 Subject: [PATCH 3/4] update test value --- test/testSource.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/testSource.cpp b/test/testSource.cpp index 6ce5cfada..4254e7f84 100644 --- a/test/testSource.cpp +++ b/test/testSource.cpp @@ -110,7 +110,7 @@ TEST(SourceSNRDistribution, simpleTest) { } R2_mean/=100000.; Z_mean/=100000.; - EXPECT_NEAR(64.4, R2_mean, 1.); + EXPECT_NEAR(100.306, R2_mean, 1); EXPECT_NEAR(0., Z_mean, 0.1); } From 3ee0447414f2e6a28d6559439a957ef23e3aadf8 Mon Sep 17 00:00:00 2001 From: Julien Date: Fri, 15 Mar 2024 11:50:34 +0100 Subject: [PATCH 4/4] add comment for alpha value --- src/Source.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Source.cpp b/src/Source.cpp index 9ff0b7f4d..f0b0e698c 100644 --- a/src/Source.cpp +++ b/src/Source.cpp @@ -416,7 +416,7 @@ double SourceSNRDistribution::fz(double z) const{ } double SourceSNRDistribution::getAlpha() const { - return alpha - 1; + return alpha - 1; // -1 to account for the R-term in the volume element dV = R * dR * dphi * dz } double SourceSNRDistribution::getBeta() const {