Skip to content
Snippets Groups Projects
Commit 15711929 authored by Carl Gwilliam's avatar Carl Gwilliam
Browse files

Allow generating particles over an annulus between two radii

parent b273fead
No related branches found
No related tags found
No related merge requests found
......@@ -40,7 +40,7 @@ def faser_pgparser():
help="Specify PDG ID of particle (note plus/minus different) or list (e.g.: --pid -13 13)")
parser.add_argument("--mass", default=105.66, type=float,
help="Specify particle mass (in MeV)")
parser.add_argument("--radius", default=100., type=float,
parser.add_argument("--radius", default=100., type=float, nargs="*",
help="Specify radius (in mm)")
parser.add_argument("--angle", default=0.005, type=float_or_none,
help="Specify angular width (in Rad)")
......@@ -55,7 +55,7 @@ def faser_pgparser():
parser.add_argument("--sampler", default="log",
help="Specify energy sampling (log, lin, const, hist, hist2D)")
parser.add_argument("--hist_name", default="log",
parser.add_argument("--hist_name", default="",
help="Specify energy sampling name for hist sampler file.root:hist")
parser.add_argument("--minE", default=10., type=float,
......
......@@ -112,6 +112,10 @@ if __name__ == '__main__':
print(f"Using pid: {args.pid} => {pidarg}")
if isinstance(args.radius, list) and len(args.radius) == 1:
args.radius = args.radius[0]
# Create the simgun dictionary
# Negative radius gives uniform sampling
# Positive radius gives Gaussian sampling
......
......@@ -35,12 +35,17 @@ class RadialPosSampler(Sampler):
def r(self):
"r position sampler"
fwhm = 2*self.radius
sig = fwhm/(2 * sqrt(2 * log(2)))
if self.radius < 0:
return sqrt(random.uniform(0, abs(self.radius**2)))
#fwhm = 2*self.radius
#sig = fwhm/(2 * sqrt(2 * log(2)))
if isinstance(self.radius, list) and len(self.radius) == 2:
# If a list of length 2, generate uniformally over an annulus from r[0] to r[1]
return sqrt(random.uniform(abs(self.radius[0]**2), abs(self.radius[1]**2)))
elif self.radius < 0:
# Generate uniformally up to |r| if r is < 0
return sqrt(random.uniform(0, abs(self.radius**2)))
else:
# Else generate as a Gaussian with widht r
x = random.gauss(0, self.radius)
y = random.gauss(0, self.radius)
return sqrt(x**2 + y**2)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment