Open bryngemark opened 1 year ago
@bryngemark @awhitbeck
If we want to mimic the implementation in TS or LYSO target, could we add both energy deposit (edep) and birds attenuation (birksFactor) in the TS hits or LYSO hits instead of adding reduced energy deposit?
As Lene Kristian said, the Birks attenuation implemented in HcalSD.cxx is rather straightforward. It calculates Birks attenuation at each step and then use it to reduce energy deposit of the step (edep) before adding the edep in the new hit.
birksFactor = 1.0 / (1.0 + birksc1_ * dedx + birksc2_ * dedx * dedx);
:
:
edep *= birksFactor;
:
:
hit.addContrib(getTrackMap().findIncident(track_id), track_id, track->GetParticleDefinition()->GetPDGEncoding(), edep, track->GetGlobalTime());
@bryngemark @awhitbeck
If trigscint is going to mimic HcalSD.cxx, codes below replace the in HcalSD.cxx's codes from LINE#44 to LINE#85 for LYSO using modified Birks attenuation.
//---------------------------------------------------------------------------------------------------
// Modified Birks' Attenuation
// ===========================
//
// We can describe the quenching effects of LYSO with
// the modified Birks' attenuation,
// using the expression and the coefficients taken from
// the paper NUCLEAR SCIENCE and TECHNOLOGY,
// Vol. 1, p.218-221 (2011) for LYSO:
//
// dL/dx = C0(dE/dx) / [1 + C1(dE/dx) + C2/(dE/dx)]
//
// with:
// C0 = 0.883
// C1 = 0.0065 gcm^-2MeV^-1
// C2 = -0.241 g^2cm^-4MeV^-2
//
// To get the "dE/dx" that appears in the formula,
// which has the dimensions
// [ dE/dx ] = MeV * cm^2 / g
// the energy deposit (in MeV) needs to be divided by the
// product of the step length (in cm) and LYSO's density
//---------------------------------------------------------------------------------------------------
G4double birksFactor(1.0);
G4double stepLength = aStep->GetStepLength() / CLHEP::cm;
// Do not apply Birks for gamma and neutron!
if (stepLength > 1.0e-6) // Check, cut if necessary.
{
G4double rho = aStep->GetPreStepPoint()->GetMaterial()->GetDensity() / (CLHEP::g / CLHEP::cm3);
G4double dedx = edep / (rho * stepLength); //[MeV*cm^2/g]
birksFactor = birksc0_/(1.0 + birksc1_*dedx + birksc2_/dedx);
if (aStep->GetTrack()->GetDefinition() == G4Gamma::GammaDefinition())
birksFactor = 1.0;
if (aStep->GetTrack()->GetDefinition() == G4Neutron::NeutronDefinition())
birksFactor = 1.0;
}
hi @joyang8caltech, welcome to github :) What would be the reason to not modify the energy deposited? The way we process simhits further down the reconstruction chain (see here) uses the deposited energy to simulate a number of PEs deposited.
If we think that the effect of Birks' law on light yield should be implemented on the light yield directly, then we probably need to think again about where exactly to plug it in. I haven't thought about this enough to be sure if we'd introduce a substantial bias to the Poisson distribution by modifying the edep rather than the PEs. I had thought about it, though, as Birks' law giving an "effective" edep which we can use to construct a Poisson for the PE yield.
If we want to keep track of the effect, I'm inclined to a) modify the edep 2) store birksFactor
as an additional member in the hit class. That way the unmodified edep can be restored if one wants to study it.
@bryngemark @awhitbeck
If trigscint is going to mimic HcalSD.cxx, codes below replace the in HcalSD.cxx's codes from LINE#44 to LINE#85 for LYSO using modified Birks attenuation. //------------------------------------------------------------------------------- ...
Please feel free to create a branch named "iss1088" to code this up in TrigScintSD.cxx :
git checkout -b iss1088
git add <files, e.g. /path/to/TrigScintSD.cxx>
git commit
<in the editor that opens, add a commit message describing your change>
git push origin iss1088
Then start a pull request as instructed, and add @awhitbeck and me and one more main developer as reviewers.
I agree, this should be implemented so that it affects the energy deposited. The number of photons is just going to be a linear function of the energy deposited, so there is no benefit that I can see to doing one over the other, and it is much easier to do the corrections where you have the information needed for the de/dx calculation.
@joyang8caltech to start out from the most updated code, check out trunk by doing
> git fetch
> git pull origin trunk
before you start developing on a branch called iss1088
Just following up some of today's TS meeting here. @joyang8caltech how do the steps that @bryngemark outlined here look to you? Do you need help getting setup with LDMX-SW? If so there are good resources here.
Hi Lauren:
Thank you for raising the issue. I’ll revisit it and get back to you.
Best regards, James
Sent from iPhone
On Sep 11, 2023, at 2:59 PM, Lauren Tompkins @.***> wrote:
Just following up some of today's TS meeting here. @joyang8caltechhttps://github.com/joyang8caltech how do the steps that @bryngemarkhttps://github.com/bryngemark outlined here look to you? Do you need help getting setup with LDMX-SW? If so there are good resources herehttps://ldmx-software.github.io/Building-ldmx-sw.html.
— Reply to this email directly, view it on GitHubhttps://github.com/LDMX-Software/ldmx-sw/issues/1088#issuecomment-1714642937, or unsubscribehttps://github.com/notifications/unsubscribe-auth/A24FEATH64MILQBZ6KVO2B3XZ6CVPANCNFSM6AAAAAAQTAULAY. You are receiving this because you were mentioned.Message ID: @.***>
Hi James,
I'm happy to chat further if it's helpful to form a game plan for this! Don't hesitate to ask.
cheers, Lauren
Lauren Tompkins
Associate Professor of Physics Stanford University Office: (+1) 650.498.0694 Cell: (+1) 630.881.7209
pronouns: she/her
From: joyang8caltech @.> Sent: Thursday, September 14, 2023 11:02 AM To: LDMX-Software/ldmx-sw @.> Cc: Lauren Tompkins @.>; Comment @.> Subject: Re: [LDMX-Software/ldmx-sw] Implement Birks' law for the trigger scintillator in framework (Issue #1088)
Hi Lauren:
Thank you for raising the issue. I’ll revisit it and get back to you.
Best regards, James
Sent from iPhone
On Sep 11, 2023, at 2:59 PM, Lauren Tompkins @.***> wrote:
Just following up some of today's TS meeting here. @joyang8caltechhttps://github.com/joyang8caltech how do the steps that @bryngemarkhttps://github.com/bryngemark outlined here look to you? Do you need help getting setup with LDMX-SW? If so there are good resources herehttps://ldmx-software.github.io/Building-ldmx-sw.html.
— Reply to this email directly, view it on GitHubhttps://github.com/LDMX-Software/ldmx-sw/issues/1088#issuecomment-1714642937, or unsubscribehttps://github.com/notifications/unsubscribe-auth/A24FEATH64MILQBZ6KVO2B3XZ6CVPANCNFSM6AAAAAAQTAULAY. You are receiving this because you were mentioned.Message ID: @.***>
— Reply to this email directly, view it on GitHubhttps://github.com/LDMX-Software/ldmx-sw/issues/1088#issuecomment-1719906504, or unsubscribehttps://github.com/notifications/unsubscribe-auth/ACYOTBSZY3PRXKGFGV4MIPLX2NBETANCNFSM6AAAAAAQTAULAY. You are receiving this because you commented.Message ID: @.***>
Is your feature request related to a problem? Please describe. We need to account for Birks' law in the TS light yield calculation, especially for LYSO. Looking at how this is implemented for HCal, it seems to be a rather straightforward addition.
Describe the solution you'd like Our implementation could be based on what's done in the HCalSD.cxx
In addition we should have a configurable flag
isLYSO_
/isPlastic_
so we know which version to apply. I also note that in the HCal implementation, Birks' law is not applied to gammas and neutrons, while I think the most correct approach would be to only apply it to charged particles?Describe alternatives you've considered You could change this inside G4 but that would mean having to have custom stuff ported whenever we change G4 version, and we can access each G4 step from within the framework, so doing it in
ldmx-sw
should be equivalent and much easier.