Realistic BTOF digitization

[ssedd1123]

Briefly, what does this PR introduce?

It's an attempt to make LGAD sensor response for Barrel TOF more realistic by doing the following,

1. Spread charge deposition from strips that are hit directly to nearby strips.
2. Perform digitization by converting charge deposited to electric pulse.
3. Digitize the pulse by converting it to TDC and ADC value. ADC propto pulse height and TDC propto time when voltage crosses certain threshold.
4. Convert the TDC and ADC value back to time and Edep by linear fit.
5. Return the BTOF hit point as "TOFBarrelRecHit". Position is estimated by either weighted sum of ADC or just center of the cell with max ADC. It's weighted average by default, but the behavior can be changed by parameters.

Noise, edge correction and time talk correction will be developed in the future.

What kind of change does this PR introduce?

• [x ] Bug fix (issue Enhanced Realism in TOF Reconstruction & Digitization #1656 )
• New feature (issue #__)
• Documentation update
• [x ] Other: Implementation of digitization and clusterization for LGADs in Barrel TOF

Please check if this PR fulfills the following:

• Tests for the changes have been added
• Documentation has been added / updated
• [x ] Changes have been communicated to collaborators

Does this PR introduce breaking changes? What changes might users need to make to their code?

Does this PR change default behavior?

It replaces the content of "TOFBarrelRecHit" with results of this new clusterization. It was originally just geant point + random smearing. I have communicated to the reconstruction group and they advise that I submit this pull request so they can study the effect of having this new digitization/clusterization algorithm. They will decide if I should save the data in a new branch or replace the origin content.

[simonge]

The issues with the functionality of the code have not been addressed yet. It might not be a problem for the majority of datasets being processed at the moment but will be a significant bug and may cause unnecessary confusion down the line for any processing of events distributed throughout time.

[ssedd1123]

The issues with the functionality of the code have not been addressed yet. It might not be a problem for the majority of datasets being processed at the moment but will be a significant bug and may cause unnecessary confusion down the line for any processing of events distributed throughout time.

I believe the main conflicts lies in the meaning of tMin. The way I intended tMin and tMax to work is that pulse generation happens at time tMin < t < tMax with no signal outside of this range. The generated pulse will be passed on to the digitization class which convert the pulse to TDC and ADC value.

You proposed that tMin = _DigitizeTime(time) where time = hit.getTime() * dd4hep::ns;, but I fail to understand why that should be the case. If tMin is indeed hit.getTime(), that means the pulse always begin at time (internal EICROC time, not global time) = 0! I don't think that should be the case. If a hit arrives at t (global time) = 5 ns (let's say), then the way I envision my class to work is that for the first 5 ns, the pulse will be flat (i.e. no voltage), and then the voltage increases and decreases according to the landau template pulse.

I understand that collision time will not be zero both due to synchronization issue between the global 100 MHz clock and the local EICROC 40 MHz clock and finite size of bunches, if that's what you want me to implement, I can do it by shifting the tMin of each event by a random number. However, the last time I communicate this idea to you, it seems like that's not what you are asking either.

Another point you raise that I cannot agree is that we cannot discard hit that arrives at t > 25 ns. To the best of my knowledge, that's really how LGAD works in reality where signal that arrives outside beyond the 40 MHz cycle is simply discarded. As a matter of fact, the pulse timing are defined by the 40 MHz clock. It's why I don't believe tMin can be related to the hit time. tMin is when a new cycle of 40 MHz clock starts regardless of when a hit arrives.

You also mentioned that I need to deal with start time alignment. That's only an issue if tMin is not a constant in an event. If tMin = 0, there is no need to worry about mismatch of bin edges because every TDC is binned with respect to the 40 MHz clock regardless of when a hit arrives. That's another reason why I don't think tMin = _DigitizeTime(time) is right.

To be fair, none of the people that I talked to at BNL can firmly confirm or correct my understanding. Since the electronics are still in the design phase, anything can change. What I've said (and programmed in this pull request) is my best guess at what should happen. I program my PulseGenerator this way because I believe it reflects how the detector works in reality.

I believe we disagree on how EICROC work, and not mistakes on programming. If we can't reach an agreement on how time works, then unfortunately there is no way I can continue the development of this digitization code. Will you be interested in a Skype/Zoom meeting sometime next week so we can resolve our misunderstanding?

[Chao1009]

Hi @ssedd1123 and @simonge, I suggest we discuss this digitization issue in the simulation WG meeting on Dec. 17 (Tuesday) at 10 am EST. Would the time work for you? https://indico.bnl.gov/event/23277/

[ssedd1123]

Hi @ssedd1123 and @simonge, I suggest we discuss this digitization issue in the simulation WG meeting on Dec. 17 (Tuesday) at 10 am EST. Would the time work for you? https://indico.bnl.gov/event/23277/

That works for me. I will prepare a presentation.

[simonge]

Hi @ssedd1123 and @simonge, I suggest we discuss this digitization issue in the simulation WG meeting on Dec. 17 (Tuesday) at 10 am EST. Would the time work for you? https://indico.bnl.gov/event/23277/

That works for me. I will prepare a presentation.

I'm happy to discuss at the meeting too, hopefully discussing in person and providing some time diagrams will clear things up. I realise you suggested this earlier too.

[Chao1009]

Thanks! The dedicated discussion has been scheduled on Dec. 17 (Tuesday) at 10 am EST: https://indico.bnl.gov/event/23277/

[veprbl]

My takeaway from the meeting is that there were the following concerns raised:

1. The 25ns window starting at t=0 works currently for single event simulations where each event contains a single bunch crossing starting at t=0, but will not work for full frame simulations. The could would need to be adapted to work for that case, and we have some old input files to test against (ref. Simon's slide 2).
2. The current code assumes that both 100 MHz EIC and 40MHz HGCROC clocks are in phase=0 at t=0. This is too optimistic, as this case occurs only once in a 5 bunch crossings. This could be addressed by storing waveforms like it's done now, but also storing strobe counters for both clocks. Those are randomly sampled in the simulation part and affect points at which ADC counts are sampled. The reconstruction would have to invert that operation.
3. The current algorithms named TOFPulseGeneration and TOFPulseDigitization contain an implementation of HGCROC/EICROC simulator. If there is a piece of code that characterizes that device (namely the part that converts currents/voltages to ADC), it should be factored out into a common implementation for re-use in other detectors.

@ssedd1123 @simonge @johnlajoie @wdconinc Please comment to correct if I'm missing anything. We also need to identify what is a minimal viable product to ship this PR and what can go into further iterations.

[ssedd1123]

My takeaway from the meeting is that there were the following concerns raised:

1. The 25ns window starting at t=0 works currently for single event simulations where each event contains a single bunch crossing starting at t=0, but will not work for full frame simulations. The could would need to be adapted to work for that case, and we have some old input files to test against (ref. Simon's slide 2).
2. The current code assumes that both 100 MHz EIC and 40MHz HGCROC clocks are in phase=0 at t=0. This is too optimistic, as this case occurs only once in a 5 bunch crossings. This could be addressed by storing waveforms like it's done now, but also storing strobe counters for both clocks. Those are randomly sampled in the simulation part and affect points at which ADC counts are sampled. The reconstruction would have to invert that operation.
3. The current algorithms named TOFPulseGeneration and TOFPulseDigitization contain an implementation of HGCROC/EICROC simulator. If there is a piece of code that characterizes that device (namely the part that converts currents/voltages to ADC), it should be factored out into a common implementation for re-use in other detectors.

@ssedd1123 @simonge @johnlajoie @wdconinc Please comment to correct if I'm missing anything. We also need to identify what is a minimal viable product to ship this PR and what can go into further iterations.

My understanding is the same as yours. One more potentially problematic issue is that the PulseGeneration class now only passes a digitized version of the puluse. It would be ideal if I can pass the original pulse and let PulseDigitization class handle the digitization. I can pass Landau parameters, but then we are restricted to using Landau pulse only.

In my opinion, the optimum way is to pass the pulse as TF1 or just a functor. The voltage will be calculated from the TF1/functor in PulseDigitization when needed. However, that would require me to pass a generic TObject. @veprbl @wdconinc Do you know if it can be done? If not, passing a digitized pulse is not the end of the world either. We can interpolate between time bins as long as the bins are not too wide.

For now, I am just going to save an extra bit of information stating which 40 MHz cycle each pulse belongs to. I will keep the pulse from TOFPulseGenerator digitized the way it is. That is the simplest and cleanest way to extend the pulse classes to record all events in the same signal frame. Once I figure out a way to pass TF1 between analysis classes, I will try to revise these pulse classes further.

[simonge]

My takeaway from the meeting is that there were the following concerns raised:

1. The 25ns window starting at t=0 works currently for single event simulations where each event contains a single bunch crossing starting at t=0, but will not work for full frame simulations. The could would need to be adapted to work for that case, and we have some old input files to test against (ref. Simon's slide 2).
2. The current code assumes that both 100 MHz EIC and 40MHz HGCROC clocks are in phase=0 at t=0. This is too optimistic, as this case occurs only once in a 5 bunch crossings. This could be addressed by storing waveforms like it's done now, but also storing strobe counters for both clocks. Those are randomly sampled in the simulation part and affect points at which ADC counts are sampled. The reconstruction would have to invert that operation.
3. The current algorithms named TOFPulseGeneration and TOFPulseDigitization contain an implementation of HGCROC/EICROC simulator. If there is a piece of code that characterizes that device (namely the part that converts currents/voltages to ADC), it should be factored out into a common implementation for re-use in other detectors.

@ssedd1123 @simonge @johnlajoie @wdconinc Please comment to correct if I'm missing anything. We also need to identify what is a minimal viable product to ship this PR and what can go into further iterations.

We got a bit carried away into the details of the multiple clocks, synchronization, bunch crossing shapes and even the information provided by the TOF detector for the physics and vertex reconstruction. Fundamentally all I want to see from a minimum viable product is something which addresses your point 1. Ideally the pulse generation wouldn't care which 40MHz clock cycle it occurred in and that everything to do with clocks gets handled by the digitization.

I'm not sure your point 2 was discussed or maybe I missed the subtlety of the point until now. Ideally this would always be handled in the MC event sample but having an option to remove potential biases from standard samples would be useful too. If a pulse isn't restricted to the 25ns window and the time is just digitized to 20ps I don't think there should be any bias based on the strobe anyway. If the 140/120ps counters are introduced in the digitization code the precision of the measurement varies over the 25ns period and maybe the deadtime too.

In my opinion, the optimum way is to pass the pulse as TF1 or just a functor. The voltage will be calculated from the TF1/functor in PulseDigitization when needed. However, that would require me to pass a generic TObject. @veprbl @wdconinc Do you know if it can be done? If not, passing a digitized pulse is not the end of the world either. We can interpolate between time bins as long as the bins are not too wide.

I think as long as the binning of the digitized pulse is equal to or smaller than the precision of the signal, passing a digitized pulse is fine, the ASIC will be sampling the signal anyway. If the pulse shape itself moves away from a purely analytical function at some point to replicate some of the random processes/noise in the detector/sensor, having a workflow with digitized pulses makes more sense to me.

[veprbl]

My takeaway from the meeting is that there were the following concerns raised:

1. The 25ns window starting at t=0 works currently for single event simulations where each event contains a single bunch crossing starting at t=0, but will not work for full frame simulations. The could would need to be adapted to work for that case, and we have some old input files to test against (ref. Simon's slide 2).
2. The current code assumes that both 100 MHz EIC and 40MHz HGCROC clocks are in phase=0 at t=0. This is too optimistic, as this case occurs only once in a 5 bunch crossings. This could be addressed by storing waveforms like it's done now, but also storing strobe counters for both clocks. Those are randomly sampled in the simulation part and affect points at which ADC counts are sampled. The reconstruction would have to invert that operation.
3. The current algorithms named TOFPulseGeneration and TOFPulseDigitization contain an implementation of HGCROC/EICROC simulator. If there is a piece of code that characterizes that device (namely the part that converts currents/voltages to ADC), it should be factored out into a common implementation for re-use in other detectors.

@ssedd1123 @simonge @johnlajoie @wdconinc Please comment to correct if I'm missing anything. We also need to identify what is a minimal viable product to ship this PR and what can go into further iterations.

My understanding is the same as yours. One more potentially problematic issue is that the PulseGeneration class now only passes a digitized version of the puluse. It would be ideal if I can pass the original pulse and let PulseDigitization class handle the digitization. I can pass Landau parameters, but then we are restricted to using Landau pulse only.

In my opinion, the optimum way is to pass the pulse as TF1 or just a functor. The voltage will be calculated from the TF1/functor in PulseDigitization when needed. However, that would require me to pass a generic TObject. @veprbl @wdconinc Do you know if it can be done? If not, passing a digitized pulse is not the end of the world either. We can interpolate between time bins as long as the bins are not too wide.

For now, I am just going to save an extra bit of information stating which 40 MHz cycle each pulse belongs to. I will keep the pulse from TOFPulseGenerator digitized the way it is. That is the simplest and cleanest way to extend the pulse classes to record all events in the same signal frame. Once I figure out a way to pass TF1 between analysis classes, I will try to revise these pulse classes further.

I was thinking about the same topic. A significant portion of the discussion is about how the signals are sampled, and I get why you would want to pass an analytical signal shape into your digitization routine. However, while the framework technically allows to pass arbitrary C++ types between factories, we chose to only pass immutable PODIO structures. The types for those have to be defined in advance and we require a careful consideration and discussion for each item that we add to our data model https://github.com/eic/EDM4eic/blob/main/edm4eic.yaml. So it would be possible to make a type for a function (e.g. represented as a string + vector of floating point parameters). However, I doubt about universality of such approach, e.g. signals may come as a result of MC simulations and not have a nice analytical expression. Would it make sense to just reuse the waveform type and pass some extra oversampled values instead?

[ssedd1123]

@veprbl @simonge @wdconinc

I believe my pull request is now ready. Here's an overview of the changes I made:

1. I renamed all "TOFPulse" classes to "LGADPulse" to make them more general.
2. Dedicated configuration classes have been introduced for BTOFChargeSharing, LGADPulseGeneration, and LGADPulseDigitization to ensure consistency.
3. The time series generated by LGADPulseGeneration now includes clock cycle information. Specifically, pulse.getTime() returns the time (in ns) at which the cycle begins. For example, pulse.getTime() == 0 for hits arriving at 0 < t < 25 ns, pulse.getTime() == 25ns for hits arriving at 25 < t < 50 ns, and so on.
4. LGADPulseDigitization has been updated to handle this additional information. Hits arriving at 0 < t < 25 ns return TDC values between 0–1023, while hits arriving at 25 < t < 50 ns return TDC values between 1024–2047, and so on.
   Both LGADPulseGeneration and LGADPulseDigitization now support handling negative time.

I believe these modifications are viable, but there are some edge cases that I should let you know. If a pulse is generated near the edges of bin cycles, and the voltage exceeds the threshold to trigger a TDC in the first EICROC cycle but peaks in the second EICROC cycle, the code will return the correct TDC value but an incorrect ADC value because the peak is being cut off. This leads to information loss. However, I don't see this as an issue because it mirrors a real-world problem that the hardware designers need to deal with. Since the electronics are still in the design phase, there is no definitive answer on how they will behave in this scenario. I propose leaving this issue open and revisiting it once we now how the electronics will behave in the real world.

Another known issue is the absence of dead time in the current implementation. Some of my colleagues are trying to find out how much dead time there will be for each cell, but, as mentioned above, the specifics are uncertain and subject to change. Dead time can be added to the implementation later once the behavior of the electronics is better understood.