Charge reconstruction in large-area photomultipliers
M. Grassi M. Montuschi M. Baldoncini F. Mantovani B. Ricci G. Andronico V. Antonelli M. Bellato E. Bernieri A. Brigatti R. Brugnera A. Budano M. Buscemi S. Bussino R. Caruso D. Chiesa D. Corti F. Dal Corso X. F. Ding S. Dusini A. Fabbri G. Fiorentini R. Ford A. Formozov G. Galet A. Garfagnini M. Giammarchi A. Giaz A. Insolia R. Isocrate I. Lippi F. Longhitano D. Lo Presti P. Lombardi F. Marini S. M. Mari C. Martellini E. Meroni M. Mezzetto L. Miramonti S. Monforte M. Nastasi F. Ortica A. Paoloni S. Parmeggiano D. Pedretti N. Pelliccia R. Pompilio E. Previtali G. Ranucci A. C. Re A. Romani P. Saggese G. Salamanna F. H. Sawy G. Settanta M. Sisti C. Sirignano M. Spinetti L. Stanco V. Strati G. Verde L. Votano
physics.ins-det













a,b,1] M. Grassi Now at APC Laboratory - IN2P3, Paris, France
c,d,2]M. MontuschiCorresponding author 
c,d]M. Baldoncini
c,d]F. Mantovani
c,d]B. Ricci


e]G. Andronico
a]V. Antonelli
f,g]M. Bellato
h,i]E. Bernieri
a,k]A.Brigatti
f,g]R. Brugnera
h,i]A. Budano
e,j]M. Buscemi
h,i]S. Bussino
e,j]R. Caruso
n,o]D. Chiesa
h,i]D. Corti
f,g]F. Dal Corso
a,p]X. F. Ding
f,g]S. Dusini
h,i]A.Fabbri
c,d]G.Fiorentini
q,a]R. Ford
a,k,l,m]A. Formozov
f]G. Galet
f,g]A. Garfagnini
a,k]M. Giammarchi
f,g]A. Giaz
e,j]A. Insolia
f,g]R. Isocrate
f,g]I. Lippi
e]F. Longhitano
e,j]D. Lo Presti
a,k]P. Lombardi
f,g]F. Marini
h,i]S. M. Mari
h,i]C. Martellini
a,k]E. Meroni
f,g]M. Mezzetto
a,k]L. Miramonti
e]S. Monforte
n,o]M. Nastasi
r,s]F. Ortica
t]A. Paoloni
a,k]S. Parmeggiano
f,g]D. Pedretti
r,s]N. Pelliccia
a,k]R. Pompilio
n,o]E. Previtali
a,k]G. Ranucci
a,k]A. C. Re
r,s]A. Romani
a,k]P. Saggese
h,i]G. Salamanna
f,g]F. H. Sawy
h,i]G. Settanta
n,o]M. Sisti
f,g]C. Sirignano
t]M. Spinetti
f,g]L. Stanco
c,u]V. Strati
e]G. Verde
t]L. Votano

[a] Istituto Nazionale di Fisica Nucleare, Sezione di Milano, 
Via Celoria 16, I-20133 Milano, Italy
[b] Institute of High Energy Physics,
 19B YuquanLu, 1000049 Beijing, P.R.China
[c] Dipartimento di Fisica e Scienze della Terra, Università di Ferrara,
 Via Saragat 1, I-44122 Ferrara, Italy
[d] Istituto Nazionale di Fisica Nucleare, Sezione di Ferrara,
 Via Saragat 1, I-44122 Ferrara, Italy

[e] Istituto Nazionale di Fisica Nucleare, Sezione di Catania,
 Via S. Sofia 62, I-95125 Catania, Italy



[f] Dipartimento di Fisica e Astronomia, Università di Padova, 
Via Marzolo 8 , I-35131 Padova, Italy
[g] Istituto Nazionale di Fisica Nucleare, Sezione di Padova, 
Via Marzolo 8 , I-35131 Padova, Italy



[h] Dipartimento di Matematica e Fisica, Università di Roma Tre,
 Via della Vasca Navale 84 , I-00146 Roma, Italy
[i] Istituto Nazionale di Fisica Nucleare, Sezione di Roma Tre, 
Via della Vasca Navale 84 , I-00146 Roma, Italy

[j] Dipartimento di Fisica e Astronomia, Università di Catania,
 Piazza Universitá 2, I-95131 Catania, Italy
[k] Dipartimento di Fisica e Astronomia, Università di Milano,
 Via Celoria 16, I-20133 Milano, Italy
[l] Joint Institute for Nuclear Research, 
141980 Dubna, Russia
[m] Lomonosov Moscow State University Skobeltsyn Institute of Nuclear Physics,
119234 Moscow, Russia


[n] Dipartimento di Fisica, Università di Milano Bicocca,
 P.zza della Scienza 3 , I-20126 Milano, Italy
[o] Istituto Nazionale di Fisica Nucleare, Sezione di Milano Bicocca, 
P.zza della Scienza 3 , I-20126 Milano, Italy

[p] Gran Sasso Science Institute, 
Via Crispi 7 , I-67100 L'Aquila, Italy
[q] SNOLAB,
 Lively, ON, P3Y 1N2 Canada

[r] Dipartimento di Chimica, Biologia e Biotecnologie, Università di Perugia, via Elce di Sotto 8 ,
 I-06123 Perugia, Italy
[s] Istituto Nazionale di Fisica Nucleare, Sezione di Perugia,
 Via Pascoli , I-06123 Perugia, Italy

[t] Istituto Nazionale di Fisica Nucleare, Laboratori Nazionali di Frascati, 
Via Fermi 40 , I-00044 Frascati (RM), 
Italy
[u] Istituto Nazionale di Fisica Nucleare, Laboratori Nazionali di Legnaro, 
Viale dell'Università 2 , I-35020 Legnaro (PD), Italy






montuschi@fe.infn.it



Large-area PhotoMultiplier Tubes () allow to efficiently instrument Liquid Scintillator (LS) neutrino detectors, where
large target masses are pivotal to compensate for neutrinos' extremely elusive nature. Depending on the detector light yield,
several scintillation photons stemming from the same neutrino interaction are 
likely to hit a single  in a few tens/hundreds of nanoseconds, resulting in several photoelectrons (s)
to pile-up at the  anode. In such scenario, the signal generated by each  is entangled
to the others, and an accurate  charge reconstruction becomes challenging. 
This manuscript describes an experimental method able to address the  charge reconstruction in the
case of large  pile-up, providing an unbiased charge estimator at the permille level up to 15 detected s.
The method is based on a signal filtering technique (Wiener filter) which suppresses the
noise due to both  and readout electronics, and on a Fourier-based deconvolution able to 
minimize the influence of signal distortions —such as an overshoot.
The analysis of simulated  waveforms shows that the slope of a linear regression modeling the relation between 
reconstructed and true charge values 
 improves from 0.769 ± 0.001 (without deconvolution) to 0.989 ± 0.001 (with deconvolution),
where unitary slope implies perfect reconstruction. A C++ implementation of the charge reconstruction algorithm is available 
online at <cit.>.





PDNet: Semantic Segmentation integrated with a Primal-Dual Network for Document binarization
    [
    
============================================================================================









§ INTRODUCTION


The consolidated Liquid Scintillator (LS) technology is driving neutrino physics into the era of precision calorimetry. The unprecedented scientific achievements of Borexino <cit.>, Daya Bay <cit.>, Double Chooz <cit.>, KamLAND <cit.> and RENO <cit.> experiments have been the trailblazer for a new generation of multi-kiloton detectors (JUNO <cit.>, Jinping <cit.>, RENO50 <cit.>, SNO+ <cit.>, ANDES <cit.>). Going beyond the open issues in particle physics, the perspective will be some spin-off in applied antineutrino physics <cit.>.

The experimental challenges of neutrino calorimetry revolve around improving both the energy and the spatial resolution of (anti)neutrino interaction detection. 
Since both resolution terms scale with the number of detected scintillation photons, maximizing the detector photocoverage —namely the sensitive area of the detector— is pivotal to improve them. In many current and future detectors, technical and budget constraints make the use of large-area  PhotoMultiplier Tubes () the only viable solution to achieve large photocoverage. In some cases, these constrains even justify a R&D program dedicated to develop a novel  technology <cit.>.  The largest  bulbs built to date are 20-inch diameter. Because of their large acceptance,
they typically detect many photoelectrons () per scintillation event, which are likely to pile up at the readout level. That is, the spacing between the s' time of arrival is lower or comparable to the width of a single-photoelectrons () pulse. 

When  pulses overlap, their identification becomes challenging. 
Especially if the pulses are affected by an overshoot (a distortion in the  output signal described in section <ref>), two subsequent s could easily be mis-reconstructed as a single one.
A biased  charge reconstruction not only compromises the linearity —and therefore the resolution— of the detector-wise energy estimator, but also threatens the time-based reconstruction of the event vertex. In a large detector, a precise knowledge of the event vertex is crucial to define a fiducial volume meant to reject background energy depositions arising from natural radioactivity. 

The relevance of this issue can be further appreciated by noting that the Daya Bay experiment, after 4 years of smooth data taking, 
developed a new readout system based on fast digitizers, and its associated waveform reconstruction, to better assess the linearity of the detector energy response <cit.>. 

The goal of this study is twofold: (i) to propose an open-source detector-independent charge reconstruction algorithm <cit.> processing the output of a generic fast digitizer (FADC) connected to a , and (ii) to define a procedure to assess the accuracy of the reconstructed charge, especially in the case of large  pile-up. While (i) is based on realistic signal and noise assumptions, and particular care was devoted to model most of the  peculiarities, (ii) is meant to allow a comparison between any charge reconstruction algorithm. It is worth mentioning that our algorithm is specifically designed to minimize those charge reconstruction biases introduced by the presence of an overshoot in the  pulses, and by the noise fluctuations embedded in the   output pulse.
 
 
 
 
 
 




§ PMT  WAVEFORM SIMULATION









We build a  waveform simulation with the aim to develop and validate the charge reconstruction algorithm with a known input signal.
To avoid to rely on the specification of a given manufacturer and/or model, we implement a parametric simulation of a generic  response. The emission of a  from the photocathode, its collection by the anode, and its subsequent amplification are simulated as an instantaneous  charge output.
In accordance to experimental data <cit.>, the charge value q  is generated randomly according to the distribution f(q) in eq. <ref>, adapted from <cit.>

    f(q) =  
     
    (1-ω) 1/σ√(2 π)exp[-(q-q_0)^2/2  σ^2] + ω/τexp- q/τ    q ≥ q_p 
    
    0     q < q_p

The weight ω determines the relative contribution of a Gaussian distribution with respect to an exponential tail, modeling the faction of under-amplified s. The width of the former (σ) is set to 30% to describe the typical uncertainty induced by the first dynode amplification (∼ 1/√(10)). The Gaussian mean value (q_0) is set to 1, as in the case of a perfect gain calibration. The exponential tail accounts for low amplitude hits due to elastically scattered and backscattered electrons from the first dynode <cit.>. The parameter q_p sets the threshold for the minimum charge to be considered in order to avoid to simulate the pedestal. The values of all the parameters used in eq. <ref> are listed in table <ref>, and the analysis of their variability is described in section <ref>. The time at which a charge output occurs is randomly sampled from a flat distribution spanning a 1 μs-long time window, and each charge-time pair is defined as a hit.



















A sequence of hits, resulting from several s impinging on the same , are shown  as instantaneous pulses in figure  <ref>A. 
The pulse position on the horizontal axis represents the hit time with respect to a reference t_0 —for example a global detector trigger—, and 
the pulse amplitude represents the  output charge.
We assume the  to be connected to a fast electronics able to sample its output at 1 GSample/s, and we describe the voltage drop
resulting from the detection of each hit
with a log-normal function (eq. <ref>) <cit.>. The AC coupling often used to split the  high-voltage from the  output signal might induce a distortion in the  output waveform. We refer to the component of the signal with a polarity opposite to the log-normal as the overshoot, and we follow <cit.> to describe it using the sum of a Gaussian (eq. <ref>) and an exponential tail (eq. <ref>). The values of the parameters used in eq. <ref>, <ref> and <ref> are listed in table <ref>. The values related to the  amplitude are chosen to obtain an approximate 0.5 mV/count resolution, resulting in a dynamic range of ∼50 when using a commercial 10-bit digitizer. The overshoot values are adapted from the Daya Bay experience <cit.>.  
The analytical shape of a -waveform 
is reported in eq. <ref> and shown in figure  <ref>B. We name this shape the   template. The signal-only waveform resulting from the hits shown is  figure  <ref>A is built by convolving each hit with the   template, and it is shown in figure  <ref>C.


    U_peak(t)=U_0 ·exp(-1/2(ln(t/τ_0)/σ_0)^2 )


    U_OS1(t)=U_1 ·exp(-1/2(t-t_1/σ_1)^2 )


    U_OS2(t)=U_2 ·1/exp(50 ns-t/10 ns) +1 ·exp(-t/τ_2)


    U(t) = U_peak + U_OS1 + U_OS2



To make the simulation more realistic, we add a noise waveform (figure  <ref>D) to the waveform built using only signal hits. The simulated noise includes a Gaussian component and some periodic components with fixed frequency. The former describes  an overall baseline offset and time-uncorrelated baseline fluctuations. The Gaussian mean (μ_N) is set to 1.5 ADC counts, and the Gaussian width (σ_N) to 1 ADC count. The fixed-frequency components describe potential noise sources embedded in, or due to, the readout circuit. We simulate 5 such components (n_FF), with a random amplitude generated flat in the 0.1-0.3 ADC counts range, and with a random frequency in the 1-480 MHz range. The ultimate waveform comprising signal hits and all noise components is shown in  figure  <ref>E.









§ CHARGE RECONSTRUCTION

























The simplest way to determine the  output charge is to sum up all the waveform samples collected in a readout window. However, this approach would embed both the noise and the overshoot into the charge measurement, resulting in a rough and biased charge estimator. The charge reconstruction algorithm described in this section and sketched in figure <ref> is meant to mitigate the role of both noise and overshoot.
It comprises two steps: (i) 
the Deconvolution Algorithm (DA), which filters
 the raw waveforms reducing the noise, and
deconvolves  templates out of the filtered waveforms; and
(ii) the Integration Algorithm (IA), which integrates the deconvolved waveforms
with the aim to determine the overall  output charge.
 
Both filtering and deconvolution are data-driven methods that need to be trained and tested. To this end, we rely on waveforms generated using the simulation package described in section <ref>. We produce two large datasets, one comprising mostly  waveforms (S dataset), and one with waveforms containing multiple hits (M dataset). In order to emulate a real experimental setup, the charge reconstruction algorithm is designed not to use the true hit information (time and charge). We unblind such parameters only when assessing its performance —as described in section <ref>. The first dataset emulates the behavior of a low-intensity LED placed in front of a  to measure its  response, where one typically gets either 0 hits or a single hit in the resulting waveform (p(0)≃0.3, p(1)≃0.6, and p(>1)≃0.1). In the second dataset, waveforms are generated such that the number of true hits follows a flat distribution in the [0-15] range. This distribution is chosen to investigate the performance of the reconstruction algorithm even in those events experiencing high pile-up.

The noise pattern present in simulated waveforms consists mostly of bin-to-bin uncorrelated fluctuations. On the contrary, a  pulse lasts for a few tens of ns. The difference between the two time scales is highlighted in figure <ref>, where the frequency spectra of a waveform with no hits (pure noise, red curve) and a waveform with 4 hits (blue curve) are compared. The former is rather flat, with the exception of the fixed-frequency noise components, which manifest as spikes in the frequency domain. The latter, on the contrary, is peaked at low frequencies, as a result of  pulses being much slower than the noise fluctuations. We exploit such a difference to suppress the noise by means of a Wiener Filter —a technique commonly employed in signal processing <cit.>. 
The filter is defined by the following kernel equation 	

    H[f]= |S[f]|^2/  |S[f]|^2 + |N[f]|^2


where S[f] and N[f] are the frequency spectra of the signal and of the noise.

In literature, this filter is referred to as the optimal linear filter for the removal of additive noise <cit.>. Namely, 
the coefficients H[f] of the Wiener Filter are calculated to minimize the average squared distance between the filter output and the desired 
signal <cit.>.
Such coefficients are then used to weight the frequency components of the waveform being processed. 


Eq. <ref> can be rearranged by dividing both the numerator and the denominator by the noise power spectrum |N[f]|^2, and by substituting the variable 
SNR[f] = |S[f]|^2 / |N[f]|^2. This manipulation yields 

    H[f] = SNR[f]/SNR[f] + 1

where SNR  is a ratio of the signal power to the noise power.

Eq. <ref> 
makes the frequency response of the filter more intuitive, being a real positive number in the range 0 ≤ H[f] ≤ 1.
Frequencies that are barely affected by noise (SNR[f] →∞) result in the filter being close to unity, hence
applying little or no attenuation to the input components. On the contrary, frequencies that are severely affected by noise 
(SNR∼ 0)  result in the filter to heavily attenuate them (H[f]∼ 0). In summary, the Wiener Filter attenuates each frequency 
component in proportion to an estimate of its signal-to-noise ratio.


The deconvolution is a technique meant to identify and resolve the presence of  pulses within a filtered waveform, with the aim to perform an unbiased  measurement of the  output charge. It is implemented as a division in the frequency domain between the waveform to be reconstructed and the template of a  pulse. 

We build the template by time-aligning and averaging 10^4  waveforms selected from the S dataset, which emulates calibration data collected by illuminating the  with a low intensity LED.

The selection determining which waveforms are to be used in the template building is based on a threshold-crossing criterion: only waveforms whose baseline-subtracted amplitude exceeds 30% of the mean  amplitude are retained. The distribution of  all the digitized samples in the S dataset (1000 samples per waveform) is shown in  figure <ref>, where the selection threshold clearly marks the transition between the noise region and the signal region. 

To further clean the   template from those noise contributions not suppressed by the averaging process, all frequencies above 120 MHz are stripped. 

The outcome of the deconvolution performed on the waveform in figure <ref>E is shown in figure <ref>. Here the deconvolved waveform is brought back to the time domain by means of an Inverse Fast Fourier transform, and it shows narrow pulses whose amplitude is proportional to the original hit charge.   However, rather than using the pulse amplitude, we find a better estimator of the hit charge to be the integral of the pulse. 

To minimize the contribution of residual noise to the pulse integral, and the possible biases arising from the ringing visible in the vicinity of the pulses (Gibbs effect), only waveform samples falling within a Region Of Interest (ROI) are summed up to yield the final charge value. A ROI is opened any time the deconvolved waveform crosses a threshold 
corresponding to 30% of the amplitude of a deconvolved  pulse. 
We define the minimal ROI to be 6 ns wide, and we extend it if at its end the waveform is still 
above the threshold. 

Following an extensive optimization procedure, we determine that our reconstruction algorithm is better suited to process waveforms with more than two hits. Below this value, the CPU power required to filter and deconvolve the waveform is not paying off in terms of charge reconstruction performance. As a consequence, we implement a waveform preselection to determine how to process each waveform. We define Δ t to be the time during which the  template stays over threshold (6 ns). When a waveform is found to be over threshold for at least 2 Δ t, it undergoes the full charge reconstruction comprising filtering and deconvolution. Otherwise, a simple integral of the bins over threshold is used to compute the reconstructed charge. 









§ RESULTS AND DISCUSSION






The performance of the charge reconstruction algorithm 
is assessed by processing a set of 10^4 simulated waveforms 
containing a random number of  in the range [0,15], namely the M  dataset. 
In particular we aim to show how the algorithm improves the precision and accuracy of the charge 
estimate with respect to a simple integration. 


To quantify the performance of the algorithm, we process 
the M dataset waveforms twice, once using the IA alone, and once using 
both the DA and the IA. For each of the two reconstruction approaches we build
a correlation plot using pairs of reconstructed and true charge values,
as shown in  figure <ref>.
An ideal reconstruction algorithm would result in the two quantities to 
be maximally correlated, and a linear regression would yield unitary slope
and null intercept.
Any deviation from such behavior is therefore to be interpreted as a bias in the reconstruction.
In particular, a non-null intercept (q) models any bias that does not depend on the number of
reconstructed , while a non-unitary slope (m)
models any -dependent bias effectively compromising the linearity of the reconstruction 
algorithm. 
To further measure how scattered the reconstructed charge values are 
with respect to the true charge values, we use the Pearson correlation coefficient (ρ).
We interpret the latter as an estimator of the charge reconstruction precision,
while we interpret the linear regression coefficients as a measurement of its accuracy.
By comparing the two plots shown in  figure <ref>, it can be noticed that the algorithm 
significantly improves both the precision and the accuracy of the charge reconstruction.
Indeed, m improves from 0.769 to 0.989, q improves from 0.540
to 0.053, and ρ improves from 0.979 to 0.988.
To make the meaning of these numbers more evident, we report that 
if waveforms with 5 and 10 true s are processed with the IA alone,
the average reconstructed charge is biased by 7% and 12% respectively.
While, in the case of DA+IA, the bias becomes negligible (at permille level)
in both cases. Figure <ref> additionally shows the distribution of the true and reconstructed charge values 
for events with a number of s ranging between 1 and 15.

We further assess the resilience of the reconstruction algorithm to possible distortions in the input waveform. 
The data taking of a typical neutrino experiment spans indeed several years —if not decades—, during 
which several aging issues might compromise the initial  performance. Such issues often affect 
both the noise level and the shape of  pulses. 
To evaluate how  the precision and the accuracy of the charge reconstruction degrades due
to variations in the  waveforms with varying shape, we produce new datasets in which all the
 input parameters are smeared by 10%, one at a time.

In analogy to the procedure described above, (i) we process all these datasets
using the full reconstruction algorithm (IA+DA); (ii) for each dataset we build the 
reco-true charge correlation plots; (iii) we use these plots to perform a linear regression and to compute the 
Pearson correlation coefficient. The discrepancy between the resulting values and the nominal values
are plotted in terms of residuals in figure <ref>. Some considerations follow.



  * Given the different role that ρ, m and q play in assessing the algorithm performance, only residuals 
within the same panel can be compared. That is, residuals here are meant to draw a hierarchy among the 
parameters describing the  shape, with the aim to show which of them affects the reconstruction algorithm the most. 


  * From the left panel it is evident that injecting s with varying amplitude (U_0) compromises
the algorithm precision. This is expected, since the  template now fails to describe the  pulses 
present in the  output waveforms.


  * The parameters defining the width of the 
 shape, namely τ_0  and σ_0, heavily affect the accuracy of the linear regression parameters.


  * The parameters shaping the overshoot (U_1, σ_1, t_1, U_2, τ_2) 
play a negligible role, as a consequence
of the successful overshoot stripping by the DA. 


  * The white noise amplitude
(σ_N) is able to introduce a non-negligible charge bias, suggesting that the implementation of
the filter has room to be improved. 




























§ CONCLUSIONS



This manuscript describes a new method to reconstruct the output charge of a  when sampled with a fast digitizer.
Its  originality stands in putting together several well established signal processing techniques with the aim to improve the charge reconstruction accuracy over a large dynamic range, in terms of both mean value of the reconstructed charge and of dispersion around the mean. 

The algorithm comprises a filtering step, a deconvolution step, and an integration step, which aim to reduce the  noise, 
to compensate for predictable distortions in the  waveform, and to infer the number of s detected by the .
Some of its features are summarized here below.




  * The filter allows to mitigate not only white noise, but also noise at fixed frequencies
 often introduced at the level of the readout electronics (e.g. due to grounding issues).



  * The deconvolution is based on a  template, which can be automatically derived from dark count 
events, effectively introducing a zero dead-time calibration procedure. In the case of an experiment 
using different  types, the algorithm can handle all of them naturally by computing a 
different template for each type. Moreover, the capability to build new templates continuously 
allows to account for any variation in the  performance over time.


  * The deconvolution allows an a-priori estimate of the time needed to reconstruct each waveform. 
The time needed to perform the most CPU-intensive operation (Fast Fourier Transform) 
depends only on the number of samples within a waveform, and not on its complexity.
On the contrary, an approach based on reconstructing the charge by means of an analytic fit 
would become slow and unreliable in case of large pile-up. 
Such consideration becomes even more relevant in detectors instrumented 
with a large number of s.



  * The deconvolution effectively reduces the undesired  features from the 
waveform (such as the overshoot), which can in principle bias the charge reconstruction, 
making the integration step much more robust.




We tested the reconstruction algorithm by analyzing waveform datasets produced with a custom-made  simulation.
The latter allowed us to simulate a generic  readout starting from an analytic  signal template.

We investigated the effect of processing waveforms with a shape different from the  template. We did it by producing
a sample of waveforms where the parameters describing  the  pulses were smeared randomly by 10% around their nominal value. The aim was to test the resilience of the algorithm to any change in the  waveform due, for instance, to aging issues.
We determined that the reconstruction algorithm is differently sensitive to the input parameters.
In particular, those describing the amplitude and the duration of a  pulse are the most likely to bias the reconstruction performance.


The overall charge reconstruction algorithm here described is geared to play a key role 
in improving the energy resolution of LS detectors with large photocoverage.
Such detectors are indeed expected to be severely affected by energy-related systematic uncertainties
stemming  from the charge reconstruction of  waveforms where several s pile up.
To provide the community with the possibility to test our reconstruction algorithm on different inputs, and to compare
its performance to different reconstruction tools, we made a C++ implementation available at <cit.>.








We are extremely grateful to S. Jetter for useful comments and suggestions. 
This work was partially supported by Ministero dell'Istruzione, dell'Università e della Ricerca (MIUR) under the scientific project PRIN 278 2012 (CPPYP7-003), by National Institute for Nuclear Physics (INFN) through the JUNO experiment and ITALian RADioactivity project (ITALRAD), by Progetto Agroalimentare Idrointelligente Aladin CUP D92I16000030009, and by Università degli Studi di Ferrara under the scientific project FIR-2017. 



99  

website
C++ implementation of the charge reconstruction algorithm described in this manuscript available at <http://www.fe.infn.it/CRA>

borexino
Borexino Collaboration, G. Alimonti et al., 



Neutrinos from the primary proton-proton fusion process in the Sun,
Nature 512 (2014) 383-386. 

<http://dx.doi.org/10.1038/nature13702> .



dayabay 
Daya Bay Collaboration, F. P. An et al., 
Observation of Electron-Antineutrino Disappearance at Daya Bay,
 108, 171803
Phys. Rev. Lett. 108  (2012) 171803. 

<http://dx.doi.org/10.1103/PhysRevLett.108.171803> .



doublechooz
Double Chooz Collaboration, Y. Abe et al.,
Improved measurements of the neutrino mixing angle θ_13  with the Double Chooz detector,
JHEP 10 (2014)  086  [Erratum ibid 02 (2015) 074 ]. 

<http://dx.doi.org/10.1007/JHEP10(2014)086> .


kamland
KamLAND Collaboration, A. Gando et al.,
Reactor on-off antineutrino measurement with KamLAND,
Phys. Rev. D 88  (2013) 033001.

 <https://doi.org/10.1103/PhysRevD.88.033001> .


reno
RENO Collaboration, J.H. Choi et al., 
Observation of Energy and Baseline Dependent Reactor Antineutrino Disappearance in the RENO Experiment,
Phys. Rev. Lett. 116 (2016)  211801. 

<https://doi.org/10.1103/PhysRevLett.116.211801> .


juno
JUNO Collaboration, Fengpeng An et al., 
Neutrino Physics with JUNO,
J. Phys. G: Nucl. Part. Phys. 43  (2016) 030401. 

<http://dx.doi.org/10.1088/0954-3899/43/3/030401> .


jinping
Wu, Yu-Cheng et al., 
Measurement of cosmic ray flux in the China JinPing underground laboratory,
 CPC 37  (2013) 086001. 

<https://doi.org/10.1088/1674-1137/37/8/086001> .


reno50
Soo-Bong Kim,
New results from RENO and prospects with RENO-50,
Nucl. Part. Phys. Proc.  265-266  (2015)  93-98 . 

<http://www.sciencedirect.com/science/article/pii/S2405601415003661?via%3Dihub> .
 

snoplus
SNO+ Collaboration, S. Andringa et al.,
Current Status and Future Prospects of the SNO+ Experiment,
High Energy Phys. 2016  (2016) 6194250.

<http://dx.doi.org/10.1155/2016/6194250> . 
 

andes
 
P. A. N. Machado, T. Muhlbeier, H. Nunokawa, and R. Zukanovich Funchal,
Potential of a neutrino detector in the ANDES underground laboratory for geophysics and astrophysics of neutrinos,
Phys. Rev. D 86 (2012) 125001. 

<https://doi.org/10.1103/PhysRevD.86.125001> .


bowden
N.S. Bowden, K.M. Heeger, P. Huber, C. Mariani, R.B. Vogelaar,  
Applied Antineutrino Physics 2015 - Conference Summary,

arXiv:1602.04759 (2016).


mcppmt
Y. Wang, 
Large Area MCP-PMT and its Application at JUNO, 
Talk at XVII  International Workshop on Neutrinos Telescopes (NeuTel 2017) , Venice, Italy, 2017, available at 
<https://agenda.infn.it/contributionDisplay.py?sessionId=15   contribId=55   confId=11857> .

dayabay_fadc
Y. Huang et al,  
The Flash ADC system and PMT waveform reconstruction for the Daya Bay Experiment,
arXiv:1707.03699 (2017). 

jetter
S. Jetter et al, 
PMT waveform modeling at the Daya Bay experiment,
CPC 36   (2012) 93-98.

<https://doi.org/10.1088/1674-1137/36/8/009> 

dossi
R. Dossi, A. Ianni, G. Ranucci, O.Ju. Smirnov,
Methods for precise photoelectron counting with photomultipliers,
Nucl. Instrum. Methods Phys. Res. A  451  (2000) 623-637

<http://dx.doi.org/10.1016/S0168-9002(00)00337-5> .


chess
J. Caravaca et al,
Experiment to demonstrate separation of Cherenkov and scintillation signals,
Phys. Rev. C 95 (2017) 055801
<https://doi.org/10.1103/PhysRevC.95.055801>


steven
Steven W. Smith , 
The Scientist and Engineer's Guide to Digital Signal Processing, 
California Technical Publishing, San Diego, U.S.A. (1998). 

wiener_description_6
S.V. Vaseghi,
Advanced Digital Signal Processing and Noise Reduction, 
John Wiley & Sons Ltd., Chichester U.K. (2000)


wiener_mit
A.V Oppenheim, G.C. Verghese, 
Signals, Systems, and Inference, 
Class notes “Introduction to Communication, Control and Signal Processing”, 
MIT 2010
<https://ocw.mit.edu/courses/electrical-engineering-and-computer-science/6-011-introduction-to-communication-control-and-signal-processing-spring-2010/readings/MIT6_011S10_chap11.pdf>




