-
Notifications
You must be signed in to change notification settings - Fork 73
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Detsim photon simulation functions #697
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Again these functions don't belong in invisible_cities/cities but in another folder. Personally my preference would be a new folder like invisible_cities/detsim or invisible_cities/optical_simulation but maybe a @mmkekic or @gonzaponte have more informed opinions.
Other than that a few clarifications on the functions and more tests needed.
def generate_s2_photons(x : np.array, | ||
el_gain : float, | ||
el_gain_sigma : float) -> np.array: | ||
"""generate number of EL-photons produced by secondary electrons that reach | ||
the EL (after drift and diffusion) | ||
""" | ||
n = len(x) | ||
nphs = np.random.normal(el_gain, el_gain_sigma, size = n) | ||
return nphs |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Should the nphs not take into account the x? What is x, actually? an array of electrons?
If you just wanted to calculate a yield value then the name should be changed and I'd just send the size rather than the full array.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
x is an array with the size of the number of secondary electrons that arrive to the EL, so any of the arrays from the output of diffuse_electrons
in PR #695. The reason why this is implemented in that way is because in the dataflow dictionary (PR #700) just the arrays appear. In my opinion this must be modified, so suggestions are welcome
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think 'x' is in this case a bad name, maybe in this case you can pass n as argument, the len of x and number_of_electrons
## same as in detsim_simulate_electrons_test | ||
@fixture(params = ["Kr83_nexus_v5_03_00_ACTIVE_7bar_3evts.MCRD.h5"]) | ||
def nexusMC_filename(request, ICDATADIR): | ||
return os.path.join(ICDATADIR, request.param) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this file is already in a session scope fixture in invisible_cities/conftest.py
, I'm not sure you gain much by redefining it
with tb.open_file(filename) as h5in: | ||
extents = pd.read_hdf(filename, 'MC/extents') | ||
event_ids = extents.evt_number | ||
hits_df = read_mchits_df(h5in, extents) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
We'll have to change this with the new readers implementation so it works for the new format too. It's in PR #693 if you want to have a look, hopefully it'll be ready soon
hits_df = read_mchits_df(h5in, extents) | ||
|
||
evt = np.random.choice(event_ids) | ||
hits = hits_df.loc[evt, :, :] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
hits_df.loc[evt, :, :]
can just be hits_df.loc[evt]
def test_pes_at_sensors_shape(random_event_hits, sensor_data): | ||
|
||
#### this could be set to other psf with detsim_get_psf functions | ||
psf = lambda x, y, z=0: x+y | ||
|
||
x_sensors, y_sensors = sensor_data["X"].values, sensor_data["Y"].values | ||
z_sensors = 100 | ||
|
||
xs, ys, zs, energies = random_event_hits | ||
photons = generate_s1_photons(energies, 10) | ||
|
||
pes = pes_at_sensors(xs, ys, zs, photons, | ||
x_sensors, y_sensors, z_sensors, | ||
psf) | ||
|
||
nsensors = len(sensor_data) | ||
nhits = len(xs) | ||
|
||
assert pes.shape == (nsensors, nhits) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Again, add individual tests for each function
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we need a test even for functions that are as simple as generate_s2_photons
, which would essentially be testing np.random.normal
?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The function is quite simple, and I guess the test can be quite simple too. Maybe it is fast to do the test anyhow.
I agree. I think this is a good time to create a simulation-focused folder where we can gather all non-generic, non-reco functions (diomira + detsim, essentially). As for the name, I don't have a strong preference, but it would like to be generic enough to include everything simulation-related. "simulation" or "detector_simulation" sound good to me, but you can come up with a better choice, probably. |
Respect the location of the code, I will go for a folder called 'detector_simulation' and include the functions there. |
85db141
to
9da4e94
Compare
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Not sure these need to be separate functions. They are one-liners that are used in exactly one place (as far as I know). Unless you feel like they add something essential I would just replace the function call by the line itself.
This functions will go in the city dataflow module. |
Following issue #691:
3. & 5. upload functions to generate scintillation (S1) and EL (S2) light and function to compute the pes at sensors based on the psf function in invisible_cities/cities/detsim_simulate_photons.py
This module contains three functions: S1 photons generation, S2 photon generation and computation of pes at sensors based on the psf.