CMS 3D CMS Logo

MaterialBudgetTree.cc
Go to the documentation of this file.
3 
5 
6 
7 MaterialBudgetTree::MaterialBudgetTree(std::shared_ptr<MaterialBudgetData> data, const std::string& filename )
8  : MaterialBudgetFormat( data )
9 {
10  theFile = std::make_unique<TFile>(filename.c_str(),"RECREATE");
11  theFile->cd();
12  book();
13 }
14 
15 
17 {
18  LogDebug("MaterialBudget") << "MaterialBudgetTree: Booking user TTree";
19  // create the TTree
20  theTree = std::make_unique<TTree>("T1","GeometryTest Tree");
21 
22  // GENERAL block
23  theTree->Branch("MB", &t_MB, "MB/F");
24  theTree->Branch("IL", &t_IL, "IL/F");
25 
26  // PARTICLE Block
27  theTree->Branch( "Particle ID", &t_ParticleID, "Particle_ID/I" );
28  theTree->Branch( "Particle Pt", &t_ParticlePt, "Particle_Pt/F" );
29  theTree->Branch( "Particle Eta", &t_ParticleEta, "Particle_Eta/F" );
30  theTree->Branch( "Particle Phi", &t_ParticlePhi, "Particle_Phi/F" );
31  theTree->Branch( "Particle Energy", &t_ParticleEnergy, "Particle_E/F" );
32  theTree->Branch( "Particle Mass", &t_ParticleMass, "Particle_M/F" );
33 
34  if( theData->allStepsON() ) {
35  theTree->Branch("Nsteps", &t_Nsteps, "Nsteps/I");
36  theTree->Branch("DeltaMB", t_DeltaMB, "DeltaMB[Nsteps]/F");
37  theTree->Branch("DeltaMB_SUP", t_DeltaMB_SUP, "DeltaMB_SUP[Nsteps]/F");
38  theTree->Branch("DeltaMB_SEN", t_DeltaMB_SEN, "DeltaMB_SEN[Nsteps]/F");
39  theTree->Branch("DeltaMB_CAB", t_DeltaMB_CAB, "DeltaMB_CAB[Nsteps]/F");
40  theTree->Branch("DeltaMB_COL", t_DeltaMB_COL, "DeltaMB_COL[Nsteps]/F");
41  theTree->Branch("DeltaMB_ELE", t_DeltaMB_ELE, "DeltaMB_ELE[Nsteps]/F");
42  theTree->Branch("DeltaMB_OTH", t_DeltaMB_OTH, "DeltaMB_OTH[Nsteps]/F");
43  theTree->Branch("DeltaMB_AIR", t_DeltaMB_AIR, "DeltaMB_AIR[Nsteps]/F");
44 
45  theTree->Branch("DeltaIL", t_DeltaIL, "DeltaIL[Nsteps]/F");
46  theTree->Branch("DeltaIL_SUP", t_DeltaIL_SUP, "DeltaIL_SUP[Nsteps]/F");
47  theTree->Branch("DeltaIL_SEN", t_DeltaIL_SEN, "DeltaIL_SEN[Nsteps]/F");
48  theTree->Branch("DeltaIL_CAB", t_DeltaIL_CAB, "DeltaIL_CAB[Nsteps]/F");
49  theTree->Branch("DeltaIL_COL", t_DeltaIL_COL, "DeltaIL_COL[Nsteps]/F");
50  theTree->Branch("DeltaIL_ELE", t_DeltaIL_ELE, "DeltaIL_ELE[Nsteps]/F");
51  theTree->Branch("DeltaIL_OTH", t_DeltaIL_OTH, "DeltaIL_OTH[Nsteps]/F");
52  theTree->Branch("DeltaIL_AIR", t_DeltaIL_AIR, "DeltaIL_AIR[Nsteps]/F");
53 
54  theTree->Branch("Initial X", t_InitialX, "Initial_X[Nsteps]/D");
55  theTree->Branch("Initial Y", t_InitialY, "Initial_Y[Nsteps]/D");
56  theTree->Branch("Initial Z", t_InitialZ, "Initial_Z[Nsteps]/D");
57 
58  theTree->Branch("Final X", t_FinalX, "Final_X[Nsteps]/D");
59  theTree->Branch("Final Y", t_FinalY, "Final_Y[Nsteps]/D");
60  theTree->Branch("Final Z", t_FinalZ, "Final_Z[Nsteps]/D");
61 
62  theTree->Branch("Volume ID", t_VolumeID, "VolumeID[Nsteps]/I");
63  theTree->Branch("Volume Name", t_VolumeName, "VolumeName[Nsteps]/C");
64  theTree->Branch("Volume Copy", t_VolumeCopy, "VolumeCopy[Nsteps]/I");
65  theTree->Branch("Volume X", t_VolumeX, "VolumeX[Nsteps]/F");
66  theTree->Branch("Volume Y", t_VolumeY, "VolumeY[Nsteps]/F");
67  theTree->Branch("Volume Z", t_VolumeZ, "VolumeZ[Nsteps]/F");
68  theTree->Branch("Volume X axis 1", t_VolumeXaxis1, "VolumeXaxis1[Nsteps]/F");
69  theTree->Branch("Volume X axis 2", t_VolumeXaxis2, "VolumeXaxis2[Nsteps]/F");
70  theTree->Branch("Volume X axis 3", t_VolumeXaxis3, "VolumeXaxis3[Nsteps]/F");
71  theTree->Branch("Volume Y axis 1", t_VolumeYaxis1, "VolumeYaxis1[Nsteps]/F");
72  theTree->Branch("Volume Y axis 2", t_VolumeYaxis2, "VolumeYaxis2[Nsteps]/F");
73  theTree->Branch("Volume Y axis 3", t_VolumeYaxis3, "VolumeYaxis3[Nsteps]/F");
74  theTree->Branch("Volume Z axis 1", t_VolumeZaxis1, "VolumeZaxis1[Nsteps]/F");
75  theTree->Branch("Volume Z axis 2", t_VolumeZaxis2, "VolumeZaxis2[Nsteps]/F");
76  theTree->Branch("Volume Z axis 3", t_VolumeZaxis3, "VolumeZaxis3[Nsteps]/F");
77 
78  theTree->Branch("Material ID", t_MaterialID, "MaterialID[Nsteps]/I");
79  theTree->Branch("Material Name", t_MaterialName, "MaterialName[Nsteps]/C");
80  theTree->Branch("Material X0", t_MaterialX0, "MaterialX0[Nsteps]/F");
81  theTree->Branch("Material Lambda0", t_MaterialLambda0, "MaterialLambda0[Nsteps]/F");
82  theTree->Branch("Material Density", t_MaterialDensity, "MaterialDensity[Nsteps]/F");
83 
84  theTree->Branch("Particle Step ID", t_ParticleStepID, "Step_ID[Nsteps]/I");
85  theTree->Branch("Particle Step Initial Pt", t_ParticleStepInitialPt, "Step_Initial_Pt[Nsteps]/F");
86  theTree->Branch("Particle Step Initial Eta", t_ParticleStepInitialEta, "Step_Initial_Eta[Nsteps]/F");
87  theTree->Branch("Particle Step Initial Phi", t_ParticleStepInitialPhi, "Step_Initial_Phi[Nsteps]/F");
88  theTree->Branch("Particle Step Initial Energy", t_ParticleStepInitialEnergy, "Step_Initial_E[Nsteps]/F");
89  theTree->Branch("Particle Step Initial Px", t_ParticleStepInitialPx, "Step_Initial_Px[Nsteps]/F");
90  theTree->Branch("Particle Step Initial Py", t_ParticleStepInitialPy, "Step_Initial_Py[Nsteps]/F");
91  theTree->Branch("Particle Step Initial Pz", t_ParticleStepInitialPz, "Step_Initial_Pz[Nsteps]/F");
92  theTree->Branch("Particle Step Initial Beta", t_ParticleStepInitialBeta, "Step_Initial_Beta[Nsteps]/F");
93  theTree->Branch("Particle Step Initial Gamma", t_ParticleStepInitialGamma, "Step_Initial_Gamma[Nsteps]/F");
94  theTree->Branch("Particle Step Initial Mass", t_ParticleStepInitialMass, "Step_Initial_Mass[Nsteps]/F");
95  theTree->Branch("Particle Step Final Pt", t_ParticleStepFinalPt, "Step_Final_Pt[Nsteps]/F");
96  theTree->Branch("Particle Step Final Eta", t_ParticleStepFinalEta, "Step_Final_Eta[Nsteps]/F");
97  theTree->Branch("Particle Step Final Phi", t_ParticleStepFinalPhi, "Step_Final_Phi[Nsteps]/F");
98  theTree->Branch("Particle Step Final Energy", t_ParticleStepFinalEnergy, "Step_Final_E[Nsteps]/F");
99  theTree->Branch("Particle Step Final Px", t_ParticleStepFinalPx, "Step_Final_Px[Nsteps]/F");
100  theTree->Branch("Particle Step Final Py", t_ParticleStepFinalPy, "Step_Final_Py[Nsteps]/F");
101  theTree->Branch("Particle Step Final Pz", t_ParticleStepFinalPz, "Step_Final_Pz[Nsteps]/F");
102  theTree->Branch("Particle Step Final Beta", t_ParticleStepFinalBeta, "Step_Final_Beta[Nsteps]/F");
103  theTree->Branch("Particle Step Final Gamma", t_ParticleStepFinalGamma, "Step_Final_Gamma[Nsteps]/F");
104  theTree->Branch("Particle Step Final Mass", t_ParticleStepFinalMass, "Step_Final_Mass[Nsteps]/F");
105  theTree->Branch("Particle Step Pre Interaction", t_ParticleStepPreInteraction, "Step_PreInteraction[Nsteps]/I");
106  theTree->Branch("Particle Step Post Interaction", t_ParticleStepPostInteraction, "Step_PostInteraction[Nsteps]/I");
107 
108  }
109 
110  LogDebug("MaterialBudget") << "MaterialBudgetTree: Booking user TTree done";
111 }
112 
113 
115 {
116 
117 }
118 
119 
121 {
122 
123 }
124 
126 {
127 
128  t_MB = theData->getTotalMB();
129  t_IL = theData->getTotalIL();
130  // t_Eta = theData->getEta();
131  // t_Phi = theData->getPhi();
132 
133  t_ParticleID = theData->getID();
134  t_ParticlePt = theData->getPt();
135  t_ParticleEta = theData->getEta();
136  t_ParticlePhi = theData->getPhi();
137  t_ParticleEnergy = theData->getEnergy();
138  t_ParticleMass = theData->getMass();
139 
140  if( theData->allStepsON() ) {
141 
142  t_Nsteps = theData->getNumberOfSteps();
143 
145 
146  edm::LogInfo("MaterialBudget") << "MaterialBudgetTree: Number of Steps into the tree " << t_Nsteps;
147 
148  for(int ii=0;ii<t_Nsteps;ii++) {
149 
150  t_DeltaMB[ii] = theData->getStepDmb(ii);
151  t_DeltaMB_SUP[ii] = theData->getSupportDmb(ii);
152  t_DeltaMB_SEN[ii] = theData->getSensitiveDmb(ii);
153  t_DeltaMB_CAB[ii] = theData->getCablesDmb(ii);
154  t_DeltaMB_COL[ii] = theData->getCoolingDmb(ii);
155  t_DeltaMB_ELE[ii] = theData->getElectronicsDmb(ii);
156  t_DeltaMB_OTH[ii] = theData->getOtherDmb(ii);
157  t_DeltaMB_AIR[ii] = theData->getAirDmb(ii);
158 
159  t_DeltaIL[ii] = theData->getStepDil(ii);
160  t_DeltaIL_SUP[ii] = theData->getSupportDil(ii);
161  t_DeltaIL_SEN[ii] = theData->getSensitiveDil(ii);
162  t_DeltaIL_CAB[ii] = theData->getCablesDil(ii);
163  t_DeltaIL_COL[ii] = theData->getCoolingDil(ii);
164  t_DeltaIL_ELE[ii] = theData->getElectronicsDil(ii);
165  t_DeltaIL_OTH[ii] = theData->getOtherDil(ii);
166  t_DeltaIL_AIR[ii] = theData->getAirDil(ii);
167 
168  t_InitialX[ii] = theData->getStepInitialX(ii);
169  t_InitialY[ii] = theData->getStepInitialY(ii);
170  t_InitialZ[ii] = theData->getStepInitialZ(ii);
171  t_FinalX[ii] = theData->getStepFinalX(ii);
172  t_FinalY[ii] = theData->getStepFinalY(ii);
173  t_FinalZ[ii] = theData->getStepFinalZ(ii);
174 
175  t_VolumeID[ii] = theData->getStepVolumeID(ii);
176  t_VolumeName[ii] = theData->getStepVolumeName(ii).c_str();
177  t_VolumeCopy[ii] = theData->getStepVolumeCopy(ii);
178  t_VolumeX[ii] = theData->getStepVolumeX(ii);
179  t_VolumeY[ii] = theData->getStepVolumeY(ii);
180  t_VolumeZ[ii] = theData->getStepVolumeZ(ii);
181  t_VolumeXaxis1[ii] = theData->getStepVolumeXaxis(ii).x();
182  t_VolumeXaxis2[ii] = theData->getStepVolumeXaxis(ii).y();
183  t_VolumeXaxis3[ii] = theData->getStepVolumeXaxis(ii).z();
184  t_VolumeYaxis1[ii] = theData->getStepVolumeYaxis(ii).x();
185  t_VolumeYaxis2[ii] = theData->getStepVolumeYaxis(ii).y();
186  t_VolumeYaxis3[ii] = theData->getStepVolumeYaxis(ii).z();
187  t_VolumeZaxis1[ii] = theData->getStepVolumeZaxis(ii).x();
188  t_VolumeZaxis2[ii] = theData->getStepVolumeZaxis(ii).y();
189  t_VolumeZaxis3[ii] = theData->getStepVolumeZaxis(ii).z();
190 
191  t_MaterialID[ii] = theData->getStepMaterialID(ii);
192  t_MaterialName[ii] = theData->getStepMaterialName(ii).c_str();
193  t_MaterialX0[ii] = theData->getStepMaterialX0(ii);
194  t_MaterialLambda0[ii] = theData->getStepMaterialLambda0(ii);
195  t_MaterialDensity[ii] = theData->getStepMaterialDensity(ii);
196 
197  t_ParticleStepID[ii] = theData->getStepID(ii);
198  t_ParticleStepInitialPt[ii] = theData->getStepInitialPt(ii);
199  t_ParticleStepInitialEta[ii] = theData->getStepInitialEta(ii);
200  t_ParticleStepInitialPhi[ii] = theData->getStepInitialPhi(ii);
201  t_ParticleStepInitialEnergy[ii] = theData->getStepInitialEnergy(ii);
202  t_ParticleStepInitialPx[ii] = theData->getStepInitialPx(ii);
203  t_ParticleStepInitialPy[ii] = theData->getStepInitialPy(ii);
204  t_ParticleStepInitialPz[ii] = theData->getStepInitialPz(ii);
205  t_ParticleStepInitialBeta[ii] = theData->getStepInitialBeta(ii);
206  t_ParticleStepInitialGamma[ii] = theData->getStepInitialGamma(ii);
207  t_ParticleStepInitialMass[ii] = theData->getStepInitialMass(ii);
208  t_ParticleStepFinalPt[ii] = theData->getStepFinalPt(ii);
209  t_ParticleStepFinalEta[ii] = theData->getStepFinalEta(ii);
210  t_ParticleStepFinalPhi[ii] = theData->getStepFinalPhi(ii);
211  t_ParticleStepFinalEnergy[ii] = theData->getStepFinalEnergy(ii);
212  t_ParticleStepFinalPx[ii] = theData->getStepFinalPx(ii);
213  t_ParticleStepFinalPy[ii] = theData->getStepFinalPy(ii);
214  t_ParticleStepFinalPz[ii] = theData->getStepFinalPz(ii);
215  t_ParticleStepFinalBeta[ii] = theData->getStepFinalBeta(ii);
216  t_ParticleStepFinalGamma[ii] = theData->getStepFinalGamma(ii);
217  t_ParticleStepFinalMass[ii] = theData->getStepFinalMass(ii);
218  t_ParticleStepPreInteraction[ii] = theData->getStepPreProcess(ii);
219  t_ParticleStepPostInteraction[ii] = theData->getStepPostProcess(ii);
220 
221  }
222  }
223 
224  theTree->Fill();
225 }
226 
227 
229 {
230 
231  // Prefered method to include any instruction
232  // once all the tracks are done
233 
234  edm::LogInfo("MaterialBudget") << "MaterialBudgetTree Writing TTree to ROOT file";
235 
236  theFile->cd();
237  theTree->Write();
238  theFile->Close();
239 }
#define LogDebug(id)
std::unique_ptr< TFile > theFile
float t_MaterialLambda0[MAXSTEPS]
float t_DeltaIL_CAB[MAXSTEPS]
float t_ParticleStepFinalEnergy[MAXSTEPS]
void fillPerStep() override
float t_VolumeXaxis3[MAXSTEPS]
float t_DeltaMB_CAB[MAXSTEPS]
float t_DeltaMB_SEN[MAXSTEPS]
float t_VolumeXaxis1[MAXSTEPS]
float t_ParticleStepFinalPt[MAXSTEPS]
MaterialBudgetTree(std::shared_ptr< MaterialBudgetData > data, const std::string &fileName)
float t_DeltaIL_OTH[MAXSTEPS]
float t_ParticleStepFinalPhi[MAXSTEPS]
float t_ParticleStepInitialPt[MAXSTEPS]
int t_MaterialID[MAXSTEPS]
float t_DeltaIL_ELE[MAXSTEPS]
float t_DeltaIL_SEN[MAXSTEPS]
float t_ParticleStepInitialPx[MAXSTEPS]
double t_InitialX[MAXSTEPS]
double t_FinalX[MAXSTEPS]
float t_DeltaMB_AIR[MAXSTEPS]
float t_VolumeYaxis3[MAXSTEPS]
float t_ParticleStepFinalPx[MAXSTEPS]
float t_ParticleStepFinalEta[MAXSTEPS]
std::unique_ptr< TTree > theTree
float t_ParticleStepInitialPhi[MAXSTEPS]
float t_VolumeZaxis1[MAXSTEPS]
float t_ParticleStepInitialPz[MAXSTEPS]
float t_ParticleStepFinalGamma[MAXSTEPS]
int t_ParticleStepPostInteraction[MAXSTEPS]
float t_DeltaMB_COL[MAXSTEPS]
int t_VolumeID[MAXSTEPS]
float t_VolumeY[MAXSTEPS]
int t_ParticleStepID[MAXSTEPS]
static const int MAXSTEPS
float t_ParticleStepFinalPy[MAXSTEPS]
float t_ParticleStepInitialPy[MAXSTEPS]
float t_DeltaMB_OTH[MAXSTEPS]
float t_ParticleStepInitialEta[MAXSTEPS]
float t_VolumeZaxis3[MAXSTEPS]
double t_InitialY[MAXSTEPS]
float t_VolumeYaxis2[MAXSTEPS]
const char * t_VolumeName[MAXSTEPS]
float t_ParticleStepInitialEnergy[MAXSTEPS]
int t_ParticleStepPreInteraction[MAXSTEPS]
float t_ParticleStepFinalMass[MAXSTEPS]
float t_ParticleStepInitialBeta[MAXSTEPS]
float t_DeltaIL_COL[MAXSTEPS]
ii
Definition: cuy.py:590
float t_MaterialDensity[MAXSTEPS]
float t_DeltaIL_AIR[MAXSTEPS]
const char * t_MaterialName[MAXSTEPS]
float t_ParticleStepFinalBeta[MAXSTEPS]
float t_VolumeZaxis2[MAXSTEPS]
float t_ParticleStepFinalPz[MAXSTEPS]
float t_ParticleStepInitialGamma[MAXSTEPS]
double t_FinalZ[MAXSTEPS]
float t_VolumeZ[MAXSTEPS]
void fillStartTrack() override
float t_DeltaMB_ELE[MAXSTEPS]
float t_DeltaMB[MAXSTEPS]
float t_VolumeX[MAXSTEPS]
void endOfRun() override
float t_DeltaMB_SUP[MAXSTEPS]
float t_DeltaIL[MAXSTEPS]
void fillEndTrack() override
char data[epos_bytes_allocation]
Definition: EPOS_Wrapper.h:82
float t_VolumeYaxis1[MAXSTEPS]
float t_DeltaIL_SUP[MAXSTEPS]
float t_ParticleStepInitialMass[MAXSTEPS]
double t_FinalY[MAXSTEPS]
float t_VolumeXaxis2[MAXSTEPS]
std::shared_ptr< MaterialBudgetData > theData
int t_VolumeCopy[MAXSTEPS]
double t_InitialZ[MAXSTEPS]
float t_MaterialX0[MAXSTEPS]