-
Notifications
You must be signed in to change notification settings - Fork 18
First ECAL implementation #120
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
Closed
Closed
Changes from all commits
Commits
Show all changes
6 commits
Select commit
Hold shift + click to select a range
54dfe81
First ECAL implementation from Nicolo Jacazio
8773086
Energy and position smearing is added
7314886
ECAL response updated, bug fixed
67dd6dd
Merge branch 'master' of github.com:AliceO2Group/DelphesO2
6a384e2
Merging with master
kharlov c8e40e1
Fix bug in ECal hit coordinate calculation
kharlov File filter
Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
There are no files selected for viewing
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -32,13 +32,14 @@ R__LOAD_LIBRARY(libDelphesO2) | |
| #include "TrackSmearer.hh" | ||
| #include "TOFLayer.hh" | ||
| #include "RICHdetector.hh" | ||
| #include "ECALdetector.hh" | ||
| #include "MIDdetector.hh" | ||
| #include "TrackUtils.hh" | ||
|
|
||
| #include "createO2tables.h" | ||
|
|
||
| // Detector parameters | ||
| const double Bz = 0.2; // [T] | ||
| const double Bz = 5.e-1; | ||
| // TOF | ||
| constexpr double tof_radius = 100.; // [cm] Radius of the TOF detector (used to compute acceptance) | ||
| const double tof_length = 200.; // [cm] Length of the TOF detector (used to compute acceptance) | ||
|
|
@@ -70,10 +71,10 @@ const double forward_rich_sigma = 7.e-3; // Resolution of the Forward RIC | |
| const char* inputFileAccMuonPID = "muonAccEffPID.root"; | ||
|
|
||
| // Simulation parameters | ||
| constexpr bool do_vertexing = true; // Vertexing with the O2 | ||
| constexpr bool enable_nuclei = true; // Nuclei LUTs | ||
| constexpr bool debug_qa = false; // Debug QA histograms | ||
| constexpr int tof_mismatch = 0; // Flag to configure the TOF mismatch running mode: 0 off, 1 create, 2 use | ||
| constexpr bool do_vertexing = true; // Vertexing with the O2 | ||
| constexpr bool enable_nuclei = true; | ||
| constexpr bool debug_qa = false; // Debug QA histograms | ||
| constexpr int tof_mismatch = 0; // Flag to configure the TOF mismatch running mode: 0 off, 1 create, 2 use | ||
|
|
||
| int createO2tables(const char* inputFile = "delphes.root", | ||
| const char* outputFile = "AODRun5.root", | ||
|
|
@@ -182,6 +183,9 @@ int createO2tables(const char* inputFile = "delphes.root", | |
| forward_rich_detector.setType(o2::delphes::RICHdetector::kForward); | ||
| forward_rich_detector.setRadiusIn(forward_rich_radius_in); | ||
|
|
||
| // ECAL detector | ||
| o2::delphes::ECALdetector ecal_detector; | ||
|
|
||
| // MID detector | ||
| o2::delphes::MIDdetector mid_detector; | ||
| const bool isMID = mid_detector.setup(inputFileAccMuonPID); | ||
|
|
@@ -198,6 +202,7 @@ int createO2tables(const char* inputFile = "delphes.root", | |
| MakeTreeO2trackExtra(); | ||
| MakeTreeO2ftof(); | ||
| MakeTreeO2rich(); | ||
| MakeTreeO2ecal(); | ||
| MakeTreeO2frich(); | ||
| MakeTreeO2mid(); | ||
| MakeTreeO2collision(); | ||
|
|
@@ -244,6 +249,7 @@ int createO2tables(const char* inputFile = "delphes.root", | |
| treeReader->ReadEntry(ientry); | ||
| constexpr float multEtaRange = 2.f; // Range in eta to count the charged particles | ||
| float dNdEta = 0.f; // Charged particle multiplicity to use in the efficiency evaluation | ||
| TLorentzVector pECAL; // 4-momentum of photon in ECAL | ||
|
|
||
| for (Int_t iparticle = 0; iparticle < particles->GetEntries(); ++iparticle) { // Loop over particles | ||
| auto particle = (GenParticle*)particles->At(iparticle); | ||
|
|
@@ -285,6 +291,24 @@ int createO2tables(const char* inputFile = "delphes.root", | |
| dNdEta += 1.f; | ||
| } | ||
| FillTree(kMcParticle); | ||
|
|
||
| // info for the ECAL | ||
| float posZ, posPhi; | ||
| if (ecal_detector.makeSignal(*particle, pECAL, posZ, posPhi)) { // to be updated 13.09.2021 | ||
|
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. The initalization of the TLorentzVector could be done once and used outside of the loop for tracks |
||
| printf("ECAL particle: pid=%d, p=(%g,%g,%g,%g), posZ=%f, posPhi=%f\n", | ||
| particle->PID, particle->Px, particle->Py, particle->Pz, particle->E,posZ,posPhi); | ||
| printf("ECAL p=(%g,%g,%g,%g)\n", | ||
| pECAL.Px(), pECAL.Py(), pECAL.Pz(), pECAL.E()); | ||
| ecal.fIndexCollisions = ientry + eventOffset; | ||
| ecal.fIndexMcParticles = TMath::Abs(iparticle + fOffsetLabel); | ||
| ecal.fPx = pECAL.Px(); | ||
| ecal.fPy = pECAL.Py(); | ||
| ecal.fPz = pECAL.Pz(); | ||
| ecal.fE = pECAL.E(); | ||
| ecal.fPosZ = posZ; | ||
| ecal.fPosPhi = posPhi; | ||
| FillTree(kA3ECAL); | ||
| } | ||
| if constexpr (debug_qa) { | ||
| if (!debugEffDenPart[particle->PID]) { | ||
| debugEffDenPart[particle->PID] = new TH1F(Form("denPart%i", particle->PID), Form("denPart%i;#it{p}_{T} (GeV/#it{c})", particle->PID), 1000, 0, 10); | ||
|
|
||
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,125 @@ | ||
| /// @author: Yuri Kharlov | ||
| /// @author: Nicolo' Jacazio <nicolo.jacazio@cern.ch> | ||
| /// @since 04/09/2021 | ||
|
|
||
| #include "ECALdetector.hh" | ||
| #include "TDatabasePDG.h" | ||
| #include "TRandom.h" | ||
| #include "TLorentzVector.h" | ||
|
|
||
| namespace o2 | ||
| { | ||
| namespace delphes | ||
| { | ||
|
|
||
| /*****************************************************************/ | ||
|
|
||
| void ECALdetector::setup(float resoEA, float resoEB, float resoEC, float resoXA, float resoXB) | ||
| { | ||
| mEnergyResolutionA = resoEA; | ||
| mEnergyResolutionB = resoEB; | ||
| mEnergyResolutionC = resoEC; | ||
| mPositionResolutionA = resoXA; | ||
| mPositionResolutionB = resoXB; | ||
| } | ||
|
|
||
| /*****************************************************************/ | ||
|
|
||
| bool ECALdetector::hasECAL(const Track& track) const | ||
| { | ||
| auto x = track.XOuter * 0.1; // [cm] | ||
| auto y = track.YOuter * 0.1; // [cm] | ||
| auto z = track.ZOuter * 0.1; // [cm] | ||
| /** check if hit **/ | ||
| bool ishit = false; | ||
| auto r = hypot(x, y); | ||
| ishit = (fabs(r - mRadius) < 0.001 && fabs(z) < mLength); | ||
| if (!ishit) | ||
| return false; | ||
| auto particle = (GenParticle*)track.Particle.GetObject(); | ||
| return true; | ||
| } | ||
|
|
||
| /*****************************************************************/ | ||
| bool ECALdetector::makeSignal(const GenParticle& particle, | ||
| TLorentzVector& p4ECAL, | ||
| float& posZ, | ||
| float& posPhi) | ||
| { | ||
| // Simulate fast response of ECAL to photons: | ||
| // take generated particle as input and calculate its smeared 4-momentum p4ECAL | ||
| // and hit coordinates posZ, posPhi | ||
|
|
||
| const int pid = particle.PID; | ||
| if (pid != 22) { // only photons are treated so far. e+- and MIPs will be added later. | ||
| return false; | ||
| } | ||
|
|
||
| TLorentzVector p4True = particle.P4(); // true 4-momentum | ||
| if (TMath::Abs(p4True.Eta()) > 4) { // ECAL acceptance is rougly limited by |eta|<4 | ||
| return false; | ||
| } | ||
|
|
||
| Float_t vtX = particle.X; // particle vertex X | ||
| Float_t vtY = particle.Y; // particle vertex Y | ||
| Float_t vtZ = particle.Z; // particle vertex Z | ||
|
|
||
| posPhi = p4True.Phi(); // azimuth angle of a photon hit | ||
| posZ = -1e6; | ||
| Double_t tanTheta = TMath::Tan(p4True.Theta()); | ||
| if (tanTheta != 0.) { | ||
| posZ = mRadius / tanTheta; // z-coodrinate of a photon hit | ||
| } | ||
|
|
||
| Double_t eSmeared = smearPhotonE(particle.E); // smeared photon energy | ||
| p4ECAL = smearPhotonP4(p4True); | ||
|
|
||
| return true; | ||
| } | ||
|
|
||
| /*****************************************************************/ | ||
| TLorentzVector ECALdetector::smearPhotonP4(const TLorentzVector pTrue) | ||
| { | ||
| // This function smears the photon 4-momentum from the true one via applying | ||
| // parametrized energy and coordinate resolution | ||
|
|
||
| // Get true energy from true 4-momentum and smear this energy | ||
| Double_t eTrue = pTrue.E(); | ||
| Double_t eSmeared = smearPhotonE(eTrue); | ||
| // Smear direction of 3-vector | ||
| Double_t phi = pTrue.Phi() + gRandom->Gaus(0., sigmaX(eTrue) / mRadius); | ||
| Double_t theta = pTrue.Theta() + gRandom->Gaus(0., sigmaX(eTrue) / mRadius); | ||
| // Calculate smeared components of 3-vector | ||
| Double_t pxSmeared = eSmeared * TMath::Cos(phi) * TMath::Sin(theta); | ||
| Double_t pySmeared = eSmeared * TMath::Sin(phi) * TMath::Sin(theta); | ||
| Double_t pzSmeared = eSmeared * TMath::Cos(theta); | ||
| // Construct new 4-momentum from smeared energy and 3-momentum | ||
| TLorentzVector pSmeared(pxSmeared, pySmeared, pzSmeared, eSmeared); | ||
| return pSmeared; | ||
| } | ||
| /*****************************************************************/ | ||
| Double_t ECALdetector::sigmaX(const Double_t eTrue) | ||
| { | ||
| // Calculate sigma of photon coordinate smearing [cm] | ||
| // E is the photon energy | ||
| Double_t dX = sqrt(mPositionResolutionA * mPositionResolutionA + mPositionResolutionB * mPositionResolutionB / eTrue); | ||
| return dX; | ||
| } | ||
| /*****************************************************************/ | ||
| Double_t ECALdetector::smearPhotonE(const Double_t eTrue) | ||
| { | ||
| // Smear a photon energy eTrue according to a Gaussian distribution with energy resolution parameters | ||
| // sigma of Gaussian smearing is calculated from parameters A,B,C and true energy | ||
|
|
||
| Double_t sigmaE = eTrue * sqrt(mEnergyResolutionA * mEnergyResolutionA / eTrue / eTrue + | ||
| mEnergyResolutionB * mEnergyResolutionB / eTrue + | ||
| mEnergyResolutionC * mEnergyResolutionC); | ||
| Double_t eSmeared = gRandom->Gaus(eTrue, sigmaE); | ||
| if (eSmeared < 0) | ||
| eSmeared = 0; | ||
| return eSmeared; | ||
| } | ||
| /*****************************************************************/ | ||
|
|
||
| } // namespace delphes | ||
| } // namespace o2 |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,43 @@ | ||
| /// @author: Roberto Preghenella | ||
| /// @email: preghenella@bo.infn.it | ||
|
|
||
| #ifndef _DelphesO2_ECALdetector_h_ | ||
| #define _DelphesO2_ECALdetector_h_ | ||
|
|
||
| #include "classes/DelphesClasses.h" | ||
|
|
||
| namespace o2 | ||
| { | ||
| namespace delphes | ||
| { | ||
|
|
||
| class ECALdetector | ||
| { | ||
|
|
||
| public: | ||
| ECALdetector() = default; | ||
| ~ECALdetector() = default; | ||
|
|
||
| void setup(float resoEA, float resoEB, float resoEC, float resoPosA, float resoPosB); | ||
| bool hasECAL(const Track& track) const; | ||
| bool makeSignal(const GenParticle& particle, TLorentzVector& pECAL, float& posZ, float& posPhi); | ||
|
|
||
| protected: | ||
| Double_t smearPhotonE(const Double_t eTrue); | ||
| Double_t sigmaX(const Double_t eTrue); | ||
| TLorentzVector smearPhotonP4(const TLorentzVector pTrue); | ||
|
|
||
| float mRadius = 120.; // ECAL barrel inner radius [cm] | ||
| float mLength = 200.; // ECAL half-length along beam axis [cm] | ||
|
|
||
| float mEnergyResolutionA = 0.002; // parameter A of energy resolution in GeV | ||
| float mEnergyResolutionB = 0.02; // parameter B of energy resolution in GeV^{1/2} | ||
| float mEnergyResolutionC = 0.10; // parameter C of energy resolution | ||
| float mPositionResolutionA = 0.15; // parameter A of coordinate resolution in cm | ||
| float mPositionResolutionB = 0.30; // parameter B of coordinate resolution in cm*GeV^{1/2} | ||
| }; | ||
|
|
||
| } // namespace delphes | ||
| } // namespace o2 | ||
|
|
||
| #endif /** _DelphesO2_ECALdetector_h_ **/ |
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
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.
This should not need to be changed