DUNE / larnd-sim

Simulation framework for a pixelated Liquid Argon TPC
Apache License 2.0
10 stars 26 forks source link

Missing channels in the light simulation #223

Closed jaafar-chakrani closed 2 months ago

jaafar-chakrani commented 2 months ago

I'm having a look at the issue of missing channels in the light simulation. It seems to come from how we're reading the LUTs. Each LUT file gives non-zero information (vis, time_dist, t0) only for 48 channels, and there is a total of 120 entries for the channels in the file (I'm not sure where this number comes from).

Here is how I checked this on NERSC:

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.backends.backend_pdf import PdfPages

# Open Module 0 LUT
lut = np.load('/global/cfs/cdirs/dune/www/data/2x2/simulation/larndsim_data/light_LUT/lightLUT_Mod0.npz')['arr']

# Look at the sum of the visibility along z-axis per channel
# I assume that LUT format is (voxel_x, voxel_y, voxel_z, op_channel_index)
# which seems to be right from larnd-sim code
pp = PdfPages('lightLUT_Mod0.pdf')
for op_channel_index in range(lut.shape[-1]):
    fig, ax = plt.subplots()
    m = ax.imshow(np.sum(lut['vis'][:, :, :, op_channel_index], axis=2))
    fig.colorbar(m, ax=ax)
    ax.set_title('channel = '+str(op_channel_index))
    plt.tight_layout()
    pp.savefig()

pp.close()

Here is the output: lightLUT_Mod0.pdf, which shows non-zero visibility only for channels 0 to 47:

I also checked in a similar way time_dist and t0 in the file, which show the same behavior.

Consequently, I'm not sure what is supposed to be in these LUTs. If we are expecting to have in each LUT file both TPCs, then something wrong happened at the generation of these LUTs. Otherwise, I would assume that the LUTs are giving the information for one TPC, in which case we need to use it for both TPCs (as is the case when using larndsim/bin/lightLUT.npz).

Assuming that the LUT gives information for one TPC only, I think there are only two places where there is an impact:

I did a larnd-sim run where I simply changed lut.shape[3] and vis_dat.shape[0] to 48 (i.e. the number of channels in one TPC, and also the number of channels that are non-zero in the LUT), and it seems to repopulate the missing light channels. As I mentioned before, with the current LUT files, the values of lut.shape[3] and vis_dat.shape[0] are equal to 120, which I'm not sure I understand where it comes from.

Again, assuming that we are only expecting one TPC in the LUT file (which appears to be the case?), the cleanest fix would be to get rid of the null entries for channels 48->120 in the LUT files. This means that lut.shape[3] and vis_dat.shape[0] will automatically be equal to 48 and would hopefully solve the issue without any change to the code. It's pretty straightforward to write a script to reprocess the LUT for this purpose.

Tagging @liviocali @YifanC @krwood @marjoleinvannuland @russellphysics

liviocali commented 2 months ago

I have just checked the original root files the LUT was written to. They have 48 channels. This means the error happened in the conversion to npz. Also 120 is the channel number in FSD looks like the script used for the conversion was set for FSD. It looks like the npz files were changed in Dec 22 which would fit with when we started to see these issues.

-rw-rw-r-- 1 yifanch dune 769M 2023-12-22 lightLUT_Mod0.npz
-rw-rw-r-- 1 yifanch dune 853M 2023-12-22 lightLUT_Mod123.npz
-rw-r--r-- 1 yifanch dune 1.1G 2023-12-22 OptSim_LUT_ArgonCube2x2_mod0_040123.root
-rw-r--r-- 1 yifanch dune 1.2G 2023-12-22 OptSim_LUT_ArgonCube2x2_mod123_050123.root

@YifanC Could you point me to the script that you used to convert them?

YifanC commented 2 months ago

@jaafar-chakrani Thank you for catching the bug and push for a solution! As you pointed out, the culprit is here. When we have (64, 128, 32, 120) with 0's in [:, :, :,48:] instead of (64, 128, 32, 48), time_profile becomes 0, and therefore there is always no signal.

  1. As this is mentioned in PR#224, the bug is caused by implicitly hardcoded channel number in a script (number of ND-LAr module light channel used instead of for 2x2). The scripts for this conversion is linked in README/Wiki. We should update this conversion script to avoid future problem.
  2. To continue use this periodic lut assumption, we should enforce a check between the lut optical channel numbers vs. number of optical channel provided by the detector configuration. We might need to be careful handling situations such as disable optical channels, missing light detectors etc. @jaafar-chakrani, I'd suggest to add this check to PR#224 already, Thank you!
  3. Currently the corresponding position of the optical channel given its id is implicitly defined by the lut, where ch0 is LCM on the top of the detector. (It goes top down from one side and then top down from the other side as the indices increase. The two sides are symmetric, so it doesn't matter if it is left first or right first.) In ndlar_flow, it goes from bottom up from one side and then bottom up from the other side as the optical channel indices increase. Either I missed a step in between, or we wrap it and unwrap it in opposite directions. I'd imagine it would be obvious using @marjoleinvannuland's event display. I'm happy to be proven wrong. Regardless, we should document how we define and use the channel id's in both larndsim (lut) and ndlar_flow.
krwood commented 2 months ago

Let's merge in the feature branch as Jaafar has it now, so we can move ahead with MiniRun5_beta2 tomorrow morning, and we can work on 2 (and maybe 3) from Yifan's list above in another branch.

I would add one more item to the list:

  1. Currently the light LUT is specified with a full path on NERSC filesystem, so the simulation is only supported on NERSC. We should consider including a catch that downloads the appropriate file from the NERSC portal for when users are running on other resources. (Of course, we certainly don't want to include the multiple GB file in the github repo itself.)
YifanC commented 2 months ago

For 4., currently if you are not runing larndsim on NERSC and don't specify custom path, it will use the default LUT in the repo. Perhaps to enable download the LUT is much more desirable implementation.

YifanC commented 2 months ago

Regarding 3., I had a sampling check based on this channel mapping with sipm_abs_pos. The sipm position checks out and matches with the LUT. I didn't understand the remapping in reco/light/mc_event_generator.py. It would be good to have some documentation explaining the logics.

YifanC commented 2 months ago

The warning with the optical channel check has been added. We still need to improve the ArgonCubeLUTSim, but given it is not within this repo, I will close this git issue for now.