FairRoot/PandaRoot
PndEmc.cxx
Go to the documentation of this file.
1 //
3 // PndEmc
4 //
5 // Filler of PndEmcPoint
6 //
7 // Created 14/08/06 by S.Spataro
8 //
10 
11 #include "PndEmc.h"
12 
13 #include "PndEmcPoint.h"
14 #include "PndEmcReader.h"
15 #include "PndEmcStructure.h"
16 #include "PndEmcGeoPar.h"
17 
18 #include "FairGeoInterface.h"
19 #include "FairGeoLoader.h"
20 #include "FairGeoRootBuilder.h"
21 #include "FairRootManager.h"
22 #include "FairVolume.h"
23 #include "FairGeoMedia.h"
24 //#include "FairGeoG3Builder.h"
25 #include "FairRuntimeDb.h"
26 #include "FairParIo.h"
27 #include "FairRun.h"
28 
29 #include "TObjArray.h"
30 #include "TClonesArray.h"
31 #include "TGeoMCGeometry.h"
32 #include "TGeoManager.h"
33 #include "TLorentzVector.h"
34 #include "TParticle.h"
35 #include "TVirtualMC.h"
36 #include "TGeoArb8.h"
37 #include "TGeoVoxelFinder.h"
38 #include "TGeoMatrix.h"
39 #include "TGeoCompositeShape.h"
40 #include "PndDetectorList.h"
41 #include "FairGeoMedium.h"
42 #include "PndStack.h"
43 #include "TString.h"
44 //#include "TGeant3.h"
45 
46 #include <iostream>
47 #include "math.h"
48 
49 using std::cout;
50 using std::endl;
51 
52 template <typename T>
53 int getsign(const T& a)
54 {
55 return (a>0 ? 1: a<0 ? -1 : 0);
56 }
57 
58 // ----- Default constructor -------------------------------------------
60  fTrackID(0),fVolumeID(0),fEventID(-1),fPos(0,0,0,0),fMom(0,0,0,0),fTime(0),fLength(0),fELoss(0),fPosIndex(0),fEmcCollection(),
61  bIsFastFsc(kFALSE), fStoreData(kTRUE), fwendcap(kFALSE), bwendcap(kFALSE), fgeoName2(""), fgeoName3(""), fgeoName4(""), MapperVersion(0)
62 {
63  fEmcCollection = new TClonesArray("PndEmcPoint");
64 
65 }
66 // -------------------------------------------------------------------------
67 
68 // ----- Standard constructor ------------------------------------------
69 
70 PndEmc::PndEmc(const char* name, Bool_t active, Bool_t fast, Bool_t storepnts):
71  FairDetector(name, active),
72  fTrackID(0),fVolumeID(0),fEventID(-1),fPos(0,0,0,0),fMom(0,0,0,0),fTime(0),fLength(0),fELoss(0),fPosIndex(0),fEmcCollection(),
73  bIsFastFsc(fast), fStoreData(storepnts), fwendcap(kFALSE), bwendcap(kFALSE), fgeoName2(""), fgeoName3(""), fgeoName4(""), MapperVersion(0)
74 {
75  fEmcCollection = new TClonesArray("PndEmcPoint");
76 
77 }
78 // -------------------------------------------------------------------------
79 
80 
81 
82 // ----- Destructor ----------------------------------------------------
84  if (fEmcCollection) {
85  fEmcCollection->Delete();
86  delete fEmcCollection;
87  }
88 
89 }
90 // -------------------------------------------------------------------------
91 
92 
93 
94 // ----- Public method Intialize ---------------------------------------
96  // Init function
97 
99  FairRun* sim = FairRun::Instance();
100  FairRuntimeDb* rtdb=sim->GetRuntimeDb();
101  PndEmcGeoPar* par=(PndEmcGeoPar*)(rtdb->getContainer("PndEmcGeoPar"));
102  par->setChanged();
103  par->setInputVersion(sim->GetRunId(),1);
104 }
105 // -------------------------------------------------------------------------
107  // Begin of the event
108 }
109 
110 
111 
112 // ----- Public method ProcessHits --------------------------------------
113 Bool_t PndEmc::ProcessHits(FairVolume* ) { // vol //[R.K.03/2017] unused variable(s)
114 
115  TString nam = gMC->CurrentVolName();
116 
117 
118 // if(gMC->IsTrackEntering()||gMC->IsTrackExiting()){
119 // printf("\n###################\n");
120 // if(gMC->IsTrackEntering()){
121 // printf("track is entering volume %s\n",nam.Data());
122 // }
123 // if(gMC->IsTrackExiting()){
124 // printf("track is exiting volume %s\n",nam.Data());
125 // }
126 // if(gMC->IsTrackInside()){
127 // printf("track is inside volume %s\n",nam.Data());
128 // }
129 // if (gMC->IsNewTrack()){
130 // printf("track is a new track %s\n", nam.Data());
131 // }
132 // TLorentzVector trackmomentum;
133 // TLorentzVector trackposition;
134 // gMC->TrackPosition(trackposition);
135 // gMC->TrackMomentum(trackmomentum);
136 // printf("Track is at x: %e, y: %e, z: %e, t: %e\n",trackposition.X(),trackposition.Y(),trackposition.Z(),trackposition.T());
137 // printf("Momentum is: px: %e, py: %e, pz: %e, E: %e\n",trackmomentum.Px(),trackmomentum.Py(),trackmomentum.Pz(),trackmomentum.E());
138 // printf("Steplength is: %e\n",gMC->TrackStep());
139 // printf("Deposited Energy: %e\n",gMC->Edep());
140 // printf("Particle Type: %i\n", gMC->TrackPid());
141 // printf("TrackID: %i\n", gMC->GetStack()->GetCurrentTrackNumber());
142 // printf("MotherID: %i\n", gMC->GetStack()->GetCurrentParentTrackNumber());
143 //
144 // // CurrentBoundaryNormal is not implemented.
145 // // Double_t x,y,z;
146 // // if(gMC->CurrentBoundaryNormal(x,y,z)){
147 // // printf("track is at boundary with normal x: %e, y: %e, z: %e\n",x,y,z);
148 // // printf("track is at an angle of %e to normal\n",trackposition.Angle(TVector3(x,y,z))*TMath::RadToDeg());
149 // // }
150 //
151 // printf("current Volume form gGeoManager is: %s\n",gGeoManager->GetCurrentVolume()->GetName());
152 // printf("step from gGeoManager is: %e\n",gGeoManager->GetStep());
153 // printf("###################\n");
154 // }
155 
156  if (gMC->Edep()<=0){
157  // skip all the points which have no energy loss (i.e. Entering)
158  // problem for MC truth!
159  // ((Idea: Check if particle was produced inside crystal or outside))
160  // ANY particle ENTERING and not being a NEW TRACK the crystal produces a hit
161  if (gMC->IsNewTrack()) return kTRUE;
162  if (!gMC->IsTrackEntering() && !gMC->IsTrackExiting()) return kTRUE;
163  }
164 
165 
166  // ---------------------------------------------------------------------------------
167  // Getting parameters for the ROOT file with geometry for Forward End-Cap.
168  // Each of the subvolume name for FwEndCap geometry in the ROOT file contains "Vol".
169  Int_t copyNoCrys, copyNoBox, copyNoSub, copyNoQuar, idCrys, idBox, idSub, idQuar, copyNo, id, nMod, nRow, nCrys;
170  copyNoCrys = copyNoBox = copyNoSub = copyNoQuar = idCrys = idBox = idSub = idQuar = copyNo = id = nMod = nRow = nCrys = -1;
171 
172  if (nam.Contains("Vol")){
173 
174  TString namCrys = gMC->CurrentVolOffName(0); // Crystal name
175  TString namBox = gMC->CurrentVolOffName(1); // Box name
176  TString namSub = gMC->CurrentVolOffName(2); // Subunit/HalfSubunit name
177  TString namQuar = gMC->CurrentVolOffName(3); // Quarter name
178 
179  if (namSub.Contains("SubunitVolFwEndCap")){//including both Subunits and HalfSubunits
180 
181  // Return the current volumeOff upward in the geometrical tree (copy number of the volume)
182  gMC->CurrentVolOffID(0,copyNoCrys);
183  gMC->CurrentVolOffID(1,copyNoBox);
184  gMC->CurrentVolOffID(2,copyNoSub);
185 
186  Int_t SubunitRow = -100; // Subunit/HalfSubunit row number
187  Int_t SubunitCol = -100; // Subunit/HalfSubunit column number
188  Int_t CrystalCol = -100; // FW EndCap crystal column number
189  Int_t CrystalRow = -100; // FW EndCap crystal row number
190  Int_t RestOfHalfSubunitRowNo[6] = {5, 6, 7, 8, 9, 9};//row number of 6 peripheral half Subunits
191  Int_t RestOfHalfSubunitColNo[6] = {8, 7, 6, 5, 3, 2};//column number of 6 peripheral half Subunits
192 
193  //determination of SubunitRow and SubunitCol:
194  if (copyNoSub >= 1 && copyNoSub <= 10){ // the middle row of 10 HalfSubunits
195  SubunitRow = 0;
196  SubunitCol = (int) pow(-1.,1+(copyNoSub-1)/5)*(4+(copyNoSub-1)%5);
197  }
198  else if (copyNoSub > 11 && copyNoSub < 24){ // the upper and lower pipe-region row of 14 HalfSubunits
199  SubunitRow = (int) pow(-1.,(copyNoSub-11)/7)*2;
200  SubunitCol = ((copyNoSub-11)%7)-3;
201  }
202  else if ((copyNoSub >= 25 && copyNoSub <= 27) || copyNoSub==57 || copyNoSub==58){ // the outermost left column of 5 HalfSubunits
203  if (copyNoSub >= 25 && copyNoSub <= 27)
204  SubunitRow = copyNoSub - 25;
205  else
206  SubunitRow = copyNoSub - 59;
207  SubunitCol = -9;
208  }
209  else if (copyNoSub >= 40 && copyNoSub <= 44){ // the outermost right column of 5 HalfSubunits
210  SubunitRow = 42 - copyNoSub;
211  SubunitCol = 9;
212  }
213  else if (copyNoSub == 28 || copyNoSub ==39 || copyNoSub == 45 || copyNoSub == 56) { // 4 single peripheral half Subunits
214  if (copyNoSub <= 39)
215  SubunitRow = RestOfHalfSubunitRowNo[0];
216  else
217  SubunitRow = -RestOfHalfSubunitRowNo[0];
218  SubunitCol = (int) pow(-1.,copyNoSub-1)*RestOfHalfSubunitColNo[0];
219  }
220  else if (copyNoSub == 29 || copyNoSub == 38 || copyNoSub == 46 || copyNoSub == 55) { // 4 single peripheral half Subunits
221  if (copyNoSub <= 38)
222  SubunitRow = RestOfHalfSubunitRowNo[1];
223  else
224  SubunitRow = -RestOfHalfSubunitRowNo[1];
225  SubunitCol = (int) pow(-1.,copyNoSub)*RestOfHalfSubunitColNo[1];
226  }
227  else if (copyNoSub == 30 || copyNoSub == 37 || copyNoSub == 47 || copyNoSub == 54) { // 4 single peripheral half Subunits
228  if (copyNoSub <= 37)
229  SubunitRow = RestOfHalfSubunitRowNo[2];
230  else
231  SubunitRow = -RestOfHalfSubunitRowNo[2];
232  SubunitCol = (int) pow(-1.,copyNoSub-1)*RestOfHalfSubunitColNo[2];
233  }
234  else if (copyNoSub == 31 || copyNoSub == 36 || copyNoSub == 48 || copyNoSub == 53) { // 4 single peripheral half Subunits
235  if (copyNoSub <= 36)
236  SubunitRow = RestOfHalfSubunitRowNo[3];
237  else
238  SubunitRow = -RestOfHalfSubunitRowNo[3];
239  SubunitCol = (int) pow(-1.,copyNoSub)*RestOfHalfSubunitColNo[3];
240  }
241  else if (copyNoSub == 32 || copyNoSub == 35 || copyNoSub == 49 || copyNoSub == 52) { // 4 single peripheral half Subunits
242  if (copyNoSub <= 35)
243  SubunitRow = RestOfHalfSubunitRowNo[4];
244  else
245  SubunitRow = -RestOfHalfSubunitRowNo[4];
246  SubunitCol = (int) pow(-1.,copyNoSub-1)*RestOfHalfSubunitColNo[4];
247  }
248  else if (copyNoSub == 33 || copyNoSub == 34 || copyNoSub == 50 || copyNoSub == 51) { // 4 single peripheral half Subunits
249  if (copyNoSub <= 34)
250  SubunitRow = RestOfHalfSubunitRowNo[5];
251  else
252  SubunitRow = -RestOfHalfSubunitRowNo[5];
253  SubunitCol = (int) pow(-1.,copyNoSub)*RestOfHalfSubunitColNo[5];
254  }
255  else if (copyNoSub >= 61 && copyNoSub <= 74){ // the middle column of full Subunits
256  SubunitRow = (int) pow(-1.,(copyNoSub-61)/7)*((copyNoSub - 61)%7 + 3);
257  SubunitCol = 0;
258  }
259  else if (copyNoSub%100 >= 1 && copyNoSub%100 <= 5){ // rest of the full Subunits
260  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
261  SubunitRow = 1;
262  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(3 + copyNoSub%100);
263  }
264  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
265  SubunitRow = -1;
266  SubunitCol = (int) pow(-1.,copyNoSub/100)*(3 + copyNoSub%100);
267  }
268  }
269  else if (copyNoSub%100 >= 6 && copyNoSub%100 <= 11){ // rest of the full Subunits
270  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
271  SubunitRow = 2;
272  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(2 + copyNoSub%100 - 5);
273  }
274  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
275  SubunitRow = -2;
276  SubunitCol = (int) pow(-1.,copyNoSub/100)*(2 + copyNoSub%100 - 5);
277  }
278  }
279  else if (copyNoSub%100 >= 12 && copyNoSub%100 <= 19){ // rest of the full Subunits
280  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
281  SubunitRow = 3;
282  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 11);
283  }
284  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
285  SubunitRow = -3;
286  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 11);
287  }
288  }
289  else if (copyNoSub%100 >= 20 && copyNoSub%100 <= 27){ // rest of the full Subunits
290  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
291  SubunitRow = 4;
292  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 19);
293  }
294  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
295  SubunitRow = -4;
296  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 19);
297  }
298  }
299  else if (copyNoSub%100 >= 28 && copyNoSub%100 <= 34){ // rest of the full Subunits
300  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
301  SubunitRow = 5;
302  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 27);
303  }
304  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
305  SubunitRow = -5;
306  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 27);
307  }
308  }
309  else if (copyNoSub%100 >= 35 && copyNoSub%100 <= 40){ // rest of the full Subunits
310  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
311  SubunitRow = 6;
312  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 34);
313  }
314  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
315  SubunitRow = -6;
316  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 34);
317  }
318  }
319  else if (copyNoSub%100 >= 41 && copyNoSub%100 <= 45){ // rest of the full Subunits
320  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
321  SubunitRow = 7;
322  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 40);
323  }
324  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
325  SubunitRow = -7;
326  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 40);
327  }
328  }
329  else if (copyNoSub%100 >= 46 && copyNoSub%100 <= 49){ // rest of the full Subunits
330  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
331  SubunitRow = 8;
332  SubunitCol = (int) pow(-1.,1+copyNoSub/100)*(copyNoSub%100 - 45);
333  }
334  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
335  SubunitRow = -8;
336  SubunitCol = (int) pow(-1.,copyNoSub/100)*(copyNoSub%100 - 45);
337  }
338  }
339  else if (copyNoSub%100 == 50 ){ // rest of the full Subunits
340  if (copyNoSub/100 == 1 || copyNoSub/100 == 4){
341  SubunitRow = 9;
342  SubunitCol = (int) pow(-1.,1+copyNoSub/100);
343  }
344  else {//here necessarily: copyNoSub/100 == 2 || copyNoSub/100 == 3
345  SubunitRow = -9;
346  SubunitCol = (int) pow(-1.,copyNoSub/100);
347  }
348  }
349 
350  Int_t flag=1; //this flag introducing is not necessary in this class but is also harmless
351  if (std::abs(SubunitRow)<=1 && std::abs(SubunitCol)<=3) flag=0; // empty copyNoSub in the beam-pipe area
352  else if (std::abs(SubunitCol)==9 && std::abs(SubunitRow)>=3) flag=0; // empty copyNoSub in the residual area
353  else if (std::abs(SubunitCol)==8 && std::abs(SubunitRow)>=6) flag=0; // -||-
354  else if (std::abs(SubunitCol)==7 && std::abs(SubunitRow)>=7) flag=0; // -||-
355  else if (std::abs(SubunitCol)==6 && std::abs(SubunitRow)>=8) flag=0; // -||-
356  else if (std::abs(SubunitCol)>=4 && std::abs(SubunitRow)==9) flag=0; // -||-
357 
358  //determination of CrystalRow and CrystalCol:
359  if (flag && copyNoSub >= 61){//determination of CrystalRow and CrystalCol for the 214 full Subunits
360  if(copyNoBox==1 || copyNoBox==4){
361  if(copyNoCrys==1 || copyNoCrys==3)
362  CrystalRow = SubunitRow*4 + (3-getsign(SubunitRow))/2;
363  else //here necessarily: copyNoCrys==2 || copyNoCrys==4
364  CrystalRow = SubunitRow*4 + (1-getsign(SubunitRow))/2;
365  }
366  else if(copyNoBox==2 || copyNoBox==3){
367  if(copyNoCrys==1 || copyNoCrys==3)
368  CrystalRow = SubunitRow*4 - (1+getsign(SubunitRow))/2;
369  else //here necessarily: copyNoCrys==2 || copyNoCrys==4
370  CrystalRow = SubunitRow*4 - (3+getsign(SubunitRow))/2;
371  }
372 
373  if (copyNoSub >= 61 && copyNoSub <= 74) {//determination of CrystalCol for the 14 middle column full Subunits
374  if(copyNoBox==1 || copyNoBox==3){
375  if(copyNoCrys==1 || copyNoCrys==4)
376  CrystalCol = 2;
377  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
378  CrystalCol = 1;
379  }
380  else if(copyNoBox==4 || copyNoBox==2){
381  if(copyNoCrys==1 || copyNoCrys==4)
382  CrystalCol = -1;
383  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
384  CrystalCol = -2;
385  }
386  }
387  else { //determination of CrystalCol for the other 4*50=200 full Subunits
388  if(copyNoBox==1 || copyNoBox==3){
389  if(copyNoCrys==1 || copyNoCrys==4)
390  CrystalCol = SubunitCol*4 + (3+getsign(SubunitCol))/2;
391  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
392  CrystalCol = SubunitCol*4 + (1+getsign(SubunitCol))/2;
393  }
394  else if(copyNoBox==4 || copyNoBox==2){
395  if(copyNoCrys==1 || copyNoCrys==4)
396  CrystalCol = SubunitCol*4 - (1-getsign(SubunitCol))/2;
397  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
398  CrystalCol = SubunitCol*4 - (3-getsign(SubunitCol))/2;
399  }
400  }
401  }
402 
403  else if (flag && copyNoSub <= 10) {//determination of CrystalRow and CrystalCol for the 10 middle-row half Subunits
404  if(copyNoCrys==1 || copyNoCrys==3)
405  CrystalRow = 1;
406  else //here necessarily: copyNoCrys==2 || copyNoCrys==4
407  CrystalRow = -1;
408 
409  if(copyNoBox==1){
410  if(copyNoCrys==1 || copyNoCrys==4)
411  CrystalCol = SubunitCol*4 + (3+getsign(SubunitCol))/2;
412  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
413  CrystalCol = SubunitCol*4 + (1+getsign(SubunitCol))/2;
414  }
415  else { //here necessarily: copyNoBox==2
416  if(copyNoCrys==1 || copyNoCrys==4)
417  CrystalCol = SubunitCol*4 - (1-getsign(SubunitCol))/2;
418  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
419  CrystalCol = SubunitCol*4 - (3-getsign(SubunitCol))/2;
420  }
421  }
422  else if (flag && copyNoSub > 11 && copyNoSub < 24) {//determination of CrystalRow and CrystalCol for the 10 upper and lower pipe-region row of half subunits
423  if(copyNoCrys==1 || copyNoCrys==3)
424  CrystalRow = (1+17*getsign(SubunitRow))/2;
425  else //here necessarily: copyNoCrys==2 || copyNoCrys==4
426  CrystalRow = (17*getsign(SubunitRow)-1)/2;
427 
428  if (SubunitCol == 0 && copyNoBox==1) {
429  if(copyNoCrys==1 || copyNoCrys==4)
430  CrystalCol = 2;
431  else
432  CrystalCol = 1;
433  }
434  else if (SubunitCol == 0) {//here necessarily: copyNoBox==2
435  if (copyNoCrys==1 || copyNoCrys==4)
436  CrystalCol = -1;
437  else
438  CrystalCol = -2;
439  }
440  else if(copyNoBox==1){
441  if(copyNoCrys==1 || copyNoCrys==4)
442  CrystalCol = SubunitCol*4 + (3+getsign(SubunitCol))/2;
443  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
444  CrystalCol = SubunitCol*4 + (1+getsign(SubunitCol))/2;
445  }
446  else { //here necessarily: copyNoBox==2
447  if(copyNoCrys==1 || copyNoCrys==4)
448  CrystalCol = SubunitCol*4 - (1-getsign(SubunitCol))/2;
449  else //here necessarily: copyNoCrys==3 || copyNoCrys==2
450  CrystalCol = SubunitCol*4 - (3-getsign(SubunitCol))/2;
451  }
452  }
453  else if (flag && ((copyNoSub >= 25 && copyNoSub <= 27) || copyNoSub == 57 || copyNoSub == 58 || (copyNoSub >= 40 && copyNoSub <= 44))) {//determination of CrystalRow and CrystalCol for the 2 outermost left and right columns of half Subunits
454  if(copyNoBox==1 && (copyNoCrys==1 || copyNoCrys==4))
455  if (copyNoSub == 25 || copyNoSub == 42)
456  CrystalRow = 2;
457  else
458  CrystalRow = SubunitRow*4 + (3+getsign(SubunitRow))/2;
459  else if (copyNoBox==1)//here necessarily: copyNoCrys==3 || copyNoCrys==2
460  if (copyNoSub == 25 || copyNoSub == 42)
461  CrystalRow = 1;
462  else
463  CrystalRow = SubunitRow*4 + (1+getsign(SubunitRow))/2;
464 
465  else if (copyNoBox==2 && (copyNoCrys==1 || copyNoCrys==4))
466  if (copyNoSub == 25 || copyNoSub == 42)
467  CrystalRow = -1;
468  else
469  CrystalRow = SubunitRow*4 - (1-getsign(SubunitRow))/2;
470  else if (copyNoBox==2)//here necessarily: copyNoCrys==3 || copyNoCrys==2 //FIXME Check Logic. use braces!
471  if (copyNoSub == 25 || copyNoSub == 42)
472  CrystalRow = -2;
473  else
474  CrystalRow = SubunitRow*4 - (3-getsign(SubunitRow))/2;
475 
476  if (copyNoCrys==1 || copyNoCrys==3)
477  if (!(copyNoSub >= 40 && copyNoSub <= 44))
478  CrystalCol = -36;
479  else
480  CrystalCol = 35;
481  else //here necessarily: copyNoCrys==4 || copyNoCrys==2
482  if (!(copyNoSub >= 40 && copyNoSub <= 44))
483  CrystalCol = -35;
484  else
485  CrystalCol = 36;
486  }
487  else if (flag && (copyNoSub == 28 || copyNoSub ==29 || copyNoSub == 38 || copyNoSub == 39 || copyNoSub == 45 || copyNoSub == 46 || copyNoSub == 55 || copyNoSub == 56)) {//determination of CrystalRow and CrystalCol for the 8 single peripheral half Subunits
488  if (copyNoBox==1 && (copyNoCrys==1 || copyNoCrys==4))
489  CrystalRow = SubunitRow*4 + (3-getsign(SubunitRow))/2;
490  else if (copyNoBox==1) //here necessarily: copyNoCrys==3 || copyNoCrys==2
491  CrystalRow = SubunitRow*4 + (1-getsign(SubunitRow))/2;
492  else if(copyNoCrys==1 || copyNoCrys==4) //here necessarily: copyNoBox==2
493  CrystalRow = SubunitRow*4 - (1+getsign(SubunitRow))/2;
494  else //here necessarily: copyNoBox==2 && (copyNoCrys==3 || copyNoCrys==2)
495  CrystalRow = SubunitRow*4 - (3+getsign(SubunitRow))/2;
496 
497  if (copyNoCrys==4 || copyNoCrys==2)
498  CrystalCol = SubunitCol*4 + (1-getsign(SubunitCol))/2;
499  else //here necessarily: copyNoCrys==3 || copyNoCrys==1
500  CrystalCol = SubunitCol*4 - (1+getsign(SubunitCol))/2;
501  }
502  else if (flag && ((copyNoSub >= 30 && copyNoSub <= 37) || (copyNoSub >=47 && copyNoSub <= 54))) {//determination of CrystalRow and CrystalCol for the 16 single peripheral half Subunits
503  if (copyNoCrys==1 || copyNoCrys==3)
504  CrystalRow = SubunitRow*4 + (1-3*getsign(SubunitRow))/2;
505  else //here necessarily: copyNoCrys==4 || copyNoCrys==2
506  CrystalRow = SubunitRow*4 - (1+3*getsign(SubunitRow))/2;
507 
508  if (copyNoBox==1 && (copyNoCrys==1 || copyNoCrys==4))
509  CrystalCol = SubunitCol*4 + (3+getsign(SubunitCol))/2;
510  else if (copyNoBox==1) //here necessarily: copyNoCrys==3 || copyNoCrys==2
511  CrystalCol = SubunitCol*4 + (1+getsign(SubunitCol))/2;
512  else if (copyNoCrys==1 || copyNoCrys==4) //here necessarily: copyNoBox==2
513  CrystalCol = SubunitCol*4 - (1-getsign(SubunitCol))/2;
514  else //here necessarily: copyNoBox==2 && (copyNoCrys==3 || copyNoCrys==2)
515  CrystalCol = SubunitCol*4 - (3-getsign(SubunitCol))/2;
516  }
517 
518  if (CrystalRow == -100 || CrystalCol == -100)
519  std::cout << "No assignment for CrystalRow and CrystalCol\n";
520 
521  nMod=3;
522  copyNo = 0;
523  nRow = CrystalRow + 37; //now nRow represents the FW EndCap crystal row number (-37<=CrystalRow<=37, CrystalRow!=0)----> 0<=nRow<=74, nRow!=37
524  nCrys = CrystalCol + 36; //now nCrys represents the FW EndCap crystal column number (-36<=CrystalCol<=36, CrystalCol!=0)----->0<=nCrys<=72, nCrys!=36
525 
526  //Text_t buffer[40];
527  //sprintf(buffer,"emc0%dr%dc%dcp%d",nMod, CrystalCol, CrystalRow, copyNoQuar);
528  //copyNo = copyNoQuar;
529  }
530 
531  /*
532  //used for the old version:
533  if (namCrys.Contains("CrystalVol_block"))
534  {
535  gMC->CurrentVolOffID(0,copyNoCrys);
536  nMod=3;
537  copyNo = copyNoCrys;
538  nRow = 999;
539  nCrys = 999;
540  }
541  */
542 
543 
544  else if (namQuar.Contains("Quarter4Vol")){
545  // ----- NEW Backward EndCap - with the FwEndCap geometry ----
546 
547  // Return the current volume off upward in the geometrical tree
548  // ID and copy number
549  idCrys = gMC->CurrentVolOffID(0,copyNoCrys);
550  idBox = gMC->CurrentVolOffID(1,copyNoBox);
551  idSub = gMC->CurrentVolOffID(2,copyNoSub)-1;
552  idQuar = gMC->CurrentVolOffID(3,copyNoQuar);
553 
554  Int_t col=0;//, k1=0; //[R.K. 01/2017] unused variable?
555  Int_t subrow=4; // 4 crystals in each subvolume
556  Int_t next=0; // starts (from the middle) next column, represents rows
557  copyNoSub-=1; // When geometry is created, copyNoSub starts from 1-13
558  // and the loop below starts from 0-12
559 
560  // Now, 26.02.2009, 3 crystals in the middle's subunit are added =>
561  // => number of Subunits for BwEncCap & straight geometry is the same
562  if((copyNoSub >= 0) && (copyNoSub <= 3)){
563  next = copyNoSub;
564  col = 0;
565  }else if((copyNoSub >= 4) && (copyNoSub <= 7)){
566  next = (copyNoSub-4);
567  col = 1;
568  }else if((copyNoSub >= 8) && (copyNoSub <= 10)){
569  next = (copyNoSub-8);
570  col = 2;
571  }else if((copyNoSub >= 11) && (copyNoSub <= 12)){
572  next = (copyNoSub-11);
573  col = 3;
574  }
575 
576  Int_t flag4=1;
577  // 26.02.2009
578  // "next" means "row", "col" means "col"
579 
580  if (next>1 && col>2) flag4=0; // corner's subunit + one below
581  if (next>2 && col>1) flag4=0; // + one to the left
582 
583  if (flag4!=0){
584  if ( (copyNoBox == 0) || (copyNoBox == 3) ){
585  if(copyNoCrys == 1 || copyNoCrys == 3){
586  nCrys = next*4 + 3;
587  }else if (copyNoCrys == 0 || copyNoCrys == 2){
588  nCrys = next*4 + 4;
589  }
590  }else if ( (copyNoBox == 1) || (copyNoBox == 2) ){
591  if(copyNoCrys == 0 || copyNoCrys == 2){
592  nCrys = next*4 + 2;
593  }else if (copyNoCrys == 1 || copyNoCrys == 3){
594  nCrys = next*4 + 1;
595  }
596  }
597  if ( (copyNoBox == 0) || (copyNoBox == 2) ){
598  if(copyNoCrys == 0 || copyNoCrys == 3){
599  nRow = subrow*col + 4;
600  }else if (copyNoCrys == 1 || copyNoCrys == 2){
601  nRow = subrow*col + 3;
602  }
603  }else if ( (copyNoBox == 1) || (copyNoBox == 3) ){
604  if(copyNoCrys == 0 || copyNoCrys == 3){
605  nRow = subrow*col + 2;
606  }else if (copyNoCrys == 1 || copyNoCrys == 2){
607  nRow = subrow*col + 1;
608  }
609  }
610  }
611  nMod=4;
612  copyNo = copyNoQuar;
613  }
614  else if (namQuar.Contains("QuarterNewVol")){
615  // ----- NEW Backward EndCap - with the FwEndCap geometry ----
616 
617  nMod=4;
618 
619  // Return the current volume off upward in the geometrical tree
620  // ID and copy number
621  idCrys = gMC->CurrentVolOffID(0,copyNoCrys);
622  idBox = gMC->CurrentVolOffID(1,copyNoBox);
623  idSub = gMC->CurrentVolOffID(2,copyNoSub)-1;
624  idQuar = gMC->CurrentVolOffID(3,copyNoQuar);
625 
626  Int_t col=0;
627  Int_t next=0; // starts (from the middle) next column, represents rows
628  copyNoSub-=1; // When geometry is created, copyNoSub starts from 1-13
629  // and the loop below starts from 0-12
630 
631  // Now, 26.02.2009, 3 crystals in the middle's subunit are added =>
632  // => number of Subunits for BwEncCap & straight geometry is the same
633  if((copyNoSub >= 0) && (copyNoSub <= 2)){
634  next = copyNoSub+1;
635  col = 0;
636  }else if((copyNoSub >= 3) && (copyNoSub <= 6)){
637  next = (copyNoSub-3);
638  col = 1;
639  }else if((copyNoSub >= 7) && (copyNoSub <= 10)){
640  next = (copyNoSub-7);
641  col = 2;
642  }else if((copyNoSub >= 11) && (copyNoSub <= 13)){
643  next = (copyNoSub-11);
644  col = 3;
645  }
646 
647  Int_t flag4=1;
648  // 26.02.2009
649  // "next" means "row", "col" means "col"
650 
651  if (next==3 && col==3) flag4=0; // corner's subunit + one below
652  if (next==0 && col==0) flag4=0; // + one to the left
653 
654  if (flag4!=0){
655  if ( (copyNoBox == 0) || (copyNoBox == 3) ){
656  if(copyNoCrys == 1 || copyNoCrys == 3){
657  nCrys = next*4 + 3;
658  }else if (copyNoCrys == 0 || copyNoCrys == 2){
659  nCrys = next*4 + 4;
660  }
661  }else if ( (copyNoBox == 1) || (copyNoBox == 2) ){
662  if(copyNoCrys == 0 || copyNoCrys == 2){
663  nCrys = next*4 + 2;
664  }else if (copyNoCrys == 1 || copyNoCrys == 3){
665  nCrys = next*4 + 1;
666  }
667  }
668  if ( (copyNoBox == 0) || (copyNoBox == 2) ){
669  if(copyNoCrys == 0 || copyNoCrys == 3){
670  nRow = col*4 + 4;
671  }else if (copyNoCrys == 1 || copyNoCrys == 2){
672  nRow = col*4 + 3;
673  }
674  }else if ( (copyNoBox == 1) || (copyNoBox == 3) ){
675  if(copyNoCrys == 0 || copyNoCrys == 3){
676  nRow = col*4 + 2;
677  }else if (copyNoCrys == 1 || copyNoCrys == 2){
678  nRow = col*4 + 1;
679  }
680  }
681  }
682  nMod=4;
683  copyNo = copyNoQuar;
684  }
685 
686  } //if (nam.Contains("Vol"))
687 
688 
689  // ----- Backward EndCap - 2017 version ----
690  if (nam.Contains("PWOCrystal") &&
691  TString(gMC->CurrentVolOffName(2)).Contains("BWECquarter") ){
692  nMod=4;
693  idCrys = gMC->CurrentVolOffID(0,nCrys);
694  idSub = gMC->CurrentVolOffID(1,nRow);
695  idQuar = gMC->CurrentVolOffID(2,copyNo);
696  //printf("nMod=%d nRow=%d copyNo=%d nCrys=%d\n",nMod,nRow,copyNo,nCrys);
697  }
698 
699 
700  // ---------------------------------------------------------------------------------
701 
702 
703 
704  fTrackID = gMC->GetStack()->GetCurrentTrackNumber(); // trk ID
705  fEventID = gMC->CurrentEvent();
706  fELoss = gMC->Edep();
707  fLength = gMC->TrackLength();
708  fTime = gMC->TrackTime();
709  gMC->TrackPosition(fPos); // cm
710  gMC->TrackMomentum(fMom); // GeV
711  //Int_t CrystalID; //[R.K.03/2017] unused variable
712  if(nam.Contains("CrystalType6")){
713  nMod=7;
714  copyNo=1;
715  TString namCrystal = gMC->CurrentVolOffName(1); //Crystal name
716  // TString namRow = gMC->CurrentVolOffName(2); // Row name
717  // namRow.Remove(0,3);
718  // if(namRow.IsDigit()){
719  // nRow=namRow.Atoi();
720  // } else {
721  // nMod=-1;
722  // }
723  //CrystalID= //[R.K.03/2017] unused variable
724  gMC->CurrentVolOffID(1,nCrys);
725  nRow =(nCrys-1)/5+1;
726  nCrys = (nCrys-1) % 5;
727  // printf("Crystal has name %s, ID %i, and copy number %i\n",namCrystal.Data(),CrystalID,nCrys);
728  if(namCrystal.Contains("CrystalType6a")){
729  nCrys=nCrys*2+1;
730  }else if(namCrystal.Contains("CrystalType6b")){
731  nCrys=nCrys*2+2;
732  }else{
733  nMod=-1;
734  }
735  }
736 
737  if (nam.BeginsWith("Crystal-")) { // added for barrel root file
738  int iAlveole, iCrystal;
739  char sign, type;
740  sscanf(nam, "Crystal-%d%c-%c%d", &iAlveole, &sign, &type, &iCrystal);
741 
742  // nMod
743  if (sign == 'p') nMod = 1;
744  else if (sign == 'm') nMod = 2;
745  else std::cout << "Error!!!" << std::endl;
746 
747  // nRow
748  nRow = (iAlveole - 1) * 4 + iCrystal;
749 
750  // copyNo
751  gMC->CurrentVolOffID(3, copyNo);
752 
753  // nCrys
754  Int_t iPhi;
755  gMC->CurrentVolOffID(0, iPhi);
756  if (nMod == 1) {
757  if (type == 'L') nCrys = (iPhi - 1) * 2 + 1;
758  if (type == 'R') nCrys = (iPhi - 1) * 2 + 2;
759  }
760  else if (nMod == 2) {
761  if (type == 'R') nCrys = (iPhi - 1) * 2 + 1;
762  if (type == 'L') nCrys = (iPhi - 1) * 2 + 2;
763  }
764  else {
765  std::cout << "Error!!!" << std::endl;
766  }
767  }
768 
769  if (nam.BeginsWith("emc"))
770  {
771 
772  sscanf(nam,"emc%dr%dc%d", &nMod, &nRow, &nCrys);
773 
774  // Crys 1-5000; copyNo 1-20; nRow 1-100, nMod 1-6
775  if ((nMod==1) || (nMod==2))
776  id = gMC->CurrentVolOffID(2,copyNo);
777  if ((nMod==3) || (nMod==4)|| (nMod==6))
778  id = gMC->CurrentVolOffID(1,copyNo);
779  // 1 -because the pad stays inside flayer4 (Emc4), so only "1" as inheritance.
780  // In barrel part one pad stays into Emc1 which stays inside Emc12 (and after Emc12 is
781  // copied and rotated -> the inheritance level is "2"
782  }
783  else if (bIsFastFsc)
784  {
785  nMod = 5;
786  nRow = abs((Int_t)(fPos.X()/11.))+1;
787  nCrys = abs((Int_t)(fPos.Y()/11.))+1;
788  if ((fPos.X()<0.) && (fPos.Y()>0.)) copyNo = 1;
789  if ((fPos.X()<0.) && (fPos.Y()<0.)) copyNo = 2;
790  if ((fPos.X()>0.) && (fPos.Y()<0.)) copyNo = 3;
791  if ((fPos.X()>0.) && (fPos.Y()>0.)) copyNo = 4;
792  }
793  else if (!bIsFastFsc )
794  {
795 
796  //Int_t ModCopy=0; //[R.K.02/2017] Unused variable
797  Int_t SupModCopy=0;
798  Int_t LocCopy=0;
799  Int_t nSupCol=0;
800  Int_t nSupRow = 0;
801  Int_t nModCol = 0;
802  Int_t nModRow = 0;
803  if (nam.Contains("FscSciVolume")){
804  gMC->CurrentVolOffID(4,SupModCopy);//upto FscSuperModuleVolume
805 
806  nSupCol = SupModCopy/100;
807  nSupRow = SupModCopy%100;
808 
809  gMC->CurrentVolOffID(2,LocCopy);//upto FscModuleVolume
810  nModCol = LocCopy%2;
811  nModRow = LocCopy/2;
812 
813  nCrys = (nSupCol - 1)*2 + nModCol + 1;
814  nRow = (nSupRow - 1)*2 + nModRow + 1;
815  nMod = 5;
816  copyNo = 1;
817 // cout<<"SupModCopy="<<SupModCopy<<endl;
818 // cout<<"LocCopy="<<LocCopy<<endl;
819 // cout<<"nRow="<<nRow<<endl;
820 // cout<<"nCrys="<<nCrys<<endl;
821 //
822 
823 // TGeant3* gMC3 = (TGeant3*) gMC;
824 // gMC3->Gpcxyz(); //a simple test
825  }
826  else if (nam.Contains("FscFiberVolume")){
827  nMod = 10; //for fibers different module
828  gMC->CurrentVolOffID(4,SupModCopy);//upto FscSuperModuleVolume
829 
830  nSupCol = SupModCopy/100;
831  nSupRow = SupModCopy%100;
832 
833  gMC->CurrentVolOffID(2,LocCopy);//upto FscModuleVolume
834  nModCol = LocCopy%2;
835  nModRow = LocCopy/2;
836 
837  nCrys = (nSupCol - 1)*2 + nModCol + 1;
838  nRow = (nSupRow - 1)*2 + nModRow + 1;
839 
840  //Int_t fiberID = //[R.K.03/2017] unused variable
841  gMC->CurrentVolOffID(1,copyNo); //copyNo - number of FIber
842 // cout<<"Fiber nRow="<<nRow<<endl;
843 // cout<<"Fiber nCrys="<<nCrys<<endl;
844 // cout<<"Fiber copyNo="<<copyNo<<endl;
845 
846  }
847  else{
848 
849  nam = gMC->CurrentVolOffName(2);
850  sscanf(nam,"emc%dr%dc%d", &nMod, &nRow, &nCrys);
851  if (nMod==5)
852  id = gMC->CurrentVolOffID(3,copyNo);
853  }
854  }
855 
856  fVolumeID = nMod*100000000 + nRow*1000000 + copyNo*10000 + nCrys;
857 
858  if (fVolumeID<0)
859  {
860  std::cout<<"Negative element index in EMC, name="<<nam<<std::endl;
861  return kTRUE;
862  }
863 
864  TVector3 pos(fPos.X(), fPos.Y(), fPos.Z());
865 
867  TVector3(fPos.X(), fPos.Y(), fPos.Z()),
868  TVector3(fMom.Px(), fMom.Py(), fMom.Pz()),
869  fTime, fLength, fELoss, nMod, nRow, nCrys, copyNo, gMC->IsTrackEntering() && !gMC->IsNewTrack(), gMC->IsTrackExiting());
870 
871  // Increment number of emc points for TParticle
872  //if (gMC->IsTrackEntering())
873  {
874  PndStack* stack = (PndStack*) gMC->GetStack();
875  stack->AddPoint(kEMC);
876  }
877 
878 // cout << "PndEmc: fTrackID= " << fTrackID
879 // << " fVolumeID= " << fVolumeID
880 // // << " fELoss= " << fELoss
881 // // << " fTime= " << fTime
882 // << " nMod= "<<nMod
883 // << " nRow= "<<nRow
884 // << " nCrys="<<nCrys
885 // << " copyNo="<<copyNo
886 // << endl;
887 
888  ResetParameters();
889 
890  return kTRUE;
891 }
892 // ----------------------------------------------------------------------------
893 
894 // ----- Public method EndOfEvent -----------------------------------------
896  if (fVerboseLevel) Print();
897  Reset();
898 }
899 // ----------------------------------------------------------------------------
900 
901 // ----- Public method Register -------------------------------------------
903 
904  FairRootManager::Instance()->Register("EmcPoint","Emc", fEmcCollection, fStoreData);
905 
906 }
907 // ----------------------------------------------------------------------------
908 
909 // ----- Public method GetCollection --------------------------------------
910 TClonesArray* PndEmc::GetCollection(Int_t iColl) const {
911  if (iColl == 0) return fEmcCollection;
912 
913  return NULL;
914 }
915 // ----------------------------------------------------------------------------
916 
917 // ----- Public method Print ----------------------------------------------
918 void PndEmc::Print() const {
919  Int_t nHits = fEmcCollection->GetEntriesFast();
920  //cout << "-I- PndEmc: " << nHits << " points registered in this event."
921  //<< endl;
922 
923  if (fVerboseLevel>1)
924  for (Int_t i=0; i<nHits; i++) (*fEmcCollection)[i]->Print();
925 }
926 // ----------------------------------------------------------------------------
927 
928 
929 
930 // ----- Public method Reset ----------------------------------------------
932  fEmcCollection->Delete();
933 
934  fPosIndex = 0;
935 }
936 // ----------------------------------------------------------------------------
937 
938 
939 // ------ Public method to enable/disable the storage of points ------------
940 
942 {
943  fStoreData=val;
944  return;
945 }
946 
947 // -------------------------------------------------------------------------
948 // guarda in FairRootManager::CopyClones
949 // ----- Public method CopyClones -----------------------------------------
950 void PndEmc::CopyClones(TClonesArray* cl1, TClonesArray* cl2, Int_t offset ) {
951  Int_t nEntries = cl1->GetEntriesFast();
952  //cout << "-I- PndEmc: " << nEntries << " entries to add." << endl;
953  TClonesArray& clref = *cl2;
954  PndEmcPoint* oldpoint = NULL;
955  for (Int_t i=0; i<nEntries; i++) {
956  oldpoint = (PndEmcPoint*) cl1->At(i);
957  Int_t index = oldpoint->GetTrackID() + offset;
958  oldpoint->SetTrackID(index);
959  new (clref[fPosIndex]) PndEmcPoint(*oldpoint);
960  fPosIndex++;
961  }
962  // cout << " -I- PndEmc: " << cl2->GetEntriesFast() << " merged entries."
963  // << endl;
964 }
965 // general function for chosing the combination of EMC geometry
966 void PndEmc::SetGeometryVersion(const Int_t GeoNumber) {
967 
968  MapperVersion =1;
969 
970  switch(GeoNumber){
971 
972  //TODO: Make emc_mechanics_and_module5_fsc.root default after resolution of overlap with Forward MDT
973  case 1: // new added for barrel EMC
974  SetGeometryFileNameQuadruple("emc_module12_2018v1.root","emc_module3_2012_new.root","emc_module4_2017.root","emc_module5_fsc.root");
975  MapperVersion = 1; // new MapperVersion for barrel EMC
976  break;
977 
978  case 2:
979  SetGeometryFileNameTriple("emc_module12.dat","emc_module3_2012_new.root","emc_module4_2017.root");
980  MapperVersion = 9;
981  break;
982 
983  case 3:
984  SetGeometryFileNameQuadruple("emc_module12.dat","emc_module3new.root","emc_module4_2017.root","emc_module5_fsc.root");
985  MapperVersion =3;
986  break;
987 
988  case 4:
989  SetGeometryFileNameTriple("emc_module12.dat","emc_module3new.root","emc_module4_2017.root");
990  MapperVersion =3;
991  break;
992 
993  case 5:
994  SetGeometryFileName("emc_module3.dat");
995  MapperVersion =2;
996  break;
997 
998  case 6:
999  SetGeometryFileName("emc_module3new.root");
1000  MapperVersion =3;
1001  break;
1002 
1003  case 7:
1004  SetGeometryFileName("emc_module3_2012_new.root");
1005  MapperVersion =4;
1006  break;
1007 
1008  case 8:
1009  SetGeometryFileName("emc_module4.dat");
1010  MapperVersion =2;
1011  break;
1012 
1013  case 9:
1014  SetGeometryFileName("emc_module4_StraightGeo26.root");
1015  MapperVersion = 9;
1016  break;
1017 
1018  case 10:
1019  SetGeometryFileName("emc_module4_StraightGeo26_Al.root");
1020  MapperVersion = 9;
1021  break;
1022 
1023  case 11:
1024  SetGeometryFileName("emc_module4_2017.root");
1025  MapperVersion = 9;
1026  break;
1027 
1028  case 112: // old (2008) version of backward emc geometry
1029  SetGeometryFileName("emc_module4_StraightGeo24.4.root");
1030  MapperVersion =8;
1031  break;
1032 
1033  case 12:
1034  SetGeometryFileName("emc_module4_StraightGeo24.4_Al2.root");
1035  MapperVersion = 9;
1036  break;
1037 
1038  case 13:
1039  SetGeometryFileName("emc_module4_FwEndCapGeo.root");
1040  MapperVersion = 9;
1041  break;
1042 
1043  case 14:
1044  SetGeometryFileName("emc_module4_FwEndCapGeo_Al.root");
1045  MapperVersion = 9;
1046  break;
1047 
1048  case 15:
1049  SetGeometryFileName("emc_module5_fsc.root");
1050  MapperVersion = 9;
1051  break;
1052 
1053  case 16:
1054  SetGeometryFileName("emc_module_5x5.dat");
1055  MapperVersion = 6;
1056  break;
1057 
1058  case 17:
1059  SetGeometryFileName("emc_proto60.root");
1060  MapperVersion = 7;
1061  break;
1062 
1063  case 18:
1064  SetGeometryFileName("emc_mechanics_and_module5_fsc.root");
1065  MapperVersion = 9;
1066  break;
1067 
1068  case 19:
1069  SetGeometryFileNameQuadruple("emc_module12.dat","emc_module3_2012_new.root","emc_module4_StraightGeo24.4.root","emc_mechanics_and_module5_fsc.root");
1070  MapperVersion = 9;
1071  break;
1072 
1073  case 20: // old version of barrel emc geometry
1074  SetGeometryFileNameQuadruple("emc_module12.dat","emc_module3_2012_new.root","emc_module4_2017.root","emc_module5_fsc.root");
1075  MapperVersion = 9;
1076  break;
1077 
1078  default:
1079  SetGeometryFileNameQuadruple("emc_module12.dat","emc_module3new.root","emc_module4_2017.root","emc_module5_fsc.root");
1080  MapperVersion = 9;
1081  break;
1082  }
1083 
1084  // store geo parameter
1085  FairRun *fRun = FairRun::Instance();
1086  FairRuntimeDb *rtdb= fRun->GetRuntimeDb();
1087  PndEmcGeoPar* par=(PndEmcGeoPar*)(rtdb->getContainer("PndEmcGeoPar"));
1088 
1090  par->SetGeometryVersion(GeoNumber);
1091  par->setChanged();
1092  //par->setInputVersion(fRun->GetRunId(),1);
1093 
1094 }
1095 
1096 void PndEmc::SetGeometryFileNameDouble(TString fname, TString fname2, Int_t fwbwchoice, TString geoVer)
1097 {
1098 
1099  SetGeometryFileName(fname, geoVer);
1100 
1101  //fgeoVer=geoVer;
1102  TString work = getenv("VMCWORKDIR");
1103 
1104  if (fwbwchoice==0) {
1105  fwendcap=kTRUE;
1106  fgeoName2=work+"/geometry/";
1107  fgeoName2+=fname2;
1108  cout <<"---> _new_ Forward End-Cap has been used: "<< fgeoName2<< endl;
1109  }else{
1110  bwendcap=kTRUE;
1111  fgeoName3=work+"/geometry/";
1112  fgeoName3+=fname2;
1113  cout <<"---> _new_ Backward End-Cap has been used: "<< fgeoName3<< endl;
1114  }
1115 }
1117 {
1118  fwendcap=kTRUE;
1119  bwendcap=kTRUE;
1120 
1121  SetGeometryFileName(fname, geoVer);
1122 
1123  //SetGeometryFileNameDouble(fname, fname2, geoVer);
1124 
1125  //fgeoVer=geoVer;
1126  TString work = getenv("VMCWORKDIR");
1127 
1128  fgeoName2=work+"/geometry/";
1129  fgeoName2+=fname2;
1130 
1131  fgeoName3=work+"/geometry/";
1132  fgeoName3+=fname3;
1133 
1134  //cout <<"fgeoName3 "<< fgeoName3<< endl;
1135 }
1136 
1137 void PndEmc::SetGeometryFileNameQuadruple(TString fname, TString fname2, TString fname3, TString fname4, TString geoVer)
1138 {
1139  fwendcap=kTRUE;
1140  bwendcap=kTRUE;
1141 
1142  SetGeometryFileName(fname, geoVer);
1143 
1144  //SetGeometryFileNameDouble(fname, fname2, geoVer);
1145 
1146  //fgeoVer=geoVer;
1147  TString work = getenv("VMCWORKDIR");
1148 
1149  fgeoName2=work+"/geometry/";
1150  fgeoName2+=fname2;
1151 
1152  fgeoName3=work+"/geometry/";
1153  fgeoName3+=fname3;
1154 
1155  fgeoName4=work+"/geometry/";
1156  fgeoName4+=fname4;
1157 
1158  //cout <<"fgeoName4 "<< fgeoName4<< endl;
1159 }
1160 
1161 
1162 // ----------------------------------------------------------------------------
1163  // ----- Public method ConstructGeometry ----------------------------------
1165  TString fileName=GetGeometryFileName();
1166 
1167  cout <<"fwendcap & bwendcap flags == "<< fwendcap <<" / "<<bwendcap<<endl;
1168 
1169  if (!fwendcap && !bwendcap){
1170 
1171  if (fileName.EndsWith(".dat")) {
1172  std::cout<< " " <<std::endl;
1173  std::cout<< " ====== EMC:: ConstructASCIIGeometry() ====== " <<std::endl;
1174  std::cout<< " ============================================= " <<std::endl;
1176  } else if(fileName.EndsWith("new.root") || fileName.EndsWith("proto60.root") || fileName.EndsWith("proto192.root")) {
1177  std::cout<< " " <<std::endl;
1178  std::cout<< " ====== EMC:: ConstructROOTGeometry() m3 === " <<std::endl;
1179  std::cout<< " ============================================ " <<std::endl;
1181  } else if(fileName.EndsWith("4_FwEndCapGeo.root") || fileName.EndsWith("4_StraightGeo26.root") || fgeoName.EndsWith("4_StraightGeo26_Al.root") || fileName.EndsWith("4_StraightGeo24.4.root") || fileName.EndsWith("4_StraightGeo24.4_Al2.root") || fileName.EndsWith("4_2017.root")) {
1182  std::cout<< " " <<std::endl;
1183  std::cout<< " ====== EMC:: ConstructROOTGeometry() m4 === " <<std::endl;
1184  std::cout<< " ============================================ " <<std::endl;
1186  } else if(fileName.EndsWith("5_fsc.root")) {
1187  std::cout<< " " <<std::endl;
1188  std::cout<< " ====== EMC:: ConstructROOTGeometry() m5 === " <<std::endl;
1189  std::cout<< " ============================================ " <<std::endl;
1191  } else {
1192  std::cout<< "Geometry format not supported " <<std::endl;
1193  }
1194  }else {
1195  if (fileName.EndsWith(".dat")) {
1196  std::cout<< " " <<std::endl;
1197  std::cout<< " ====== EMC 2):: ConstructASCIIGeometry() === " <<std::endl;
1198  std::cout<< " ============================================= " <<std::endl;
1200  } else if(fileName.EndsWith("2018v1.root")) {
1201  std::cout<< " " <<std::endl;
1202  std::cout<< " ====== EMC 2):: ConstructROOTGeometry() === " <<std::endl;
1203  std::cout<< " ============================================= " <<std::endl;
1205  }
1206  else {
1207  std::cout<< "You do not provide an ASCII file " <<std::endl;
1208  }
1209 
1210  Bool_t bEmc3=kFALSE, bEmc4=kFALSE, bEmc5=kFALSE;
1211  if (fgeoName2.EndsWith("new.root")) {
1212  std::cout<< " " <<std::endl;
1213  std::cout<< " ====== EMC:: ConstructRootGeometry() m3a === " <<std::endl;
1214  std::cout<< " ============================================= " <<std::endl;
1216  bEmc3 = kTRUE;
1217  }
1218  if(fgeoName3.EndsWith("4_FwEndCapGeo.root") || fgeoName3.EndsWith("4_StraightGeo26.root") || fgeoName3.EndsWith("4_StraightGeo26_Al.root") || fgeoName3.EndsWith("4_StraightGeo24.4.root") || fgeoName3.EndsWith("4_StraightGeo24.4_Al2.root") || fgeoName3.EndsWith("4_2017.root") ) {
1219  std::cout<< " " <<std::endl;
1220  std::cout<< " ====== EMC:: ConstructRootGeometry() m4a === " <<std::endl;
1221  std::cout<< " ============================================= " <<std::endl;
1222  std::cout<< "fgeoName3:: "<<fgeoName3 <<std::endl;
1224  bEmc4 = kTRUE;
1225  }
1226  if(fgeoName4.EndsWith("5_fsc.root")) {
1227  std::cout<< " " <<std::endl;
1228  std::cout<< " ====== EMC:: ConstructRootGeometry() m5a === " <<std::endl;
1229  std::cout<< " ============================================= " <<std::endl;
1230  std::cout<< "fgeoName4:: "<<fgeoName4 <<std::endl;
1232  bEmc5 = kTRUE;
1233  }
1234  if(!bEmc3 && !bEmc4 && !bEmc5 ) {
1235  std::cout<< "You do not provide a ROOT file " <<std::endl;
1236  }
1237  }
1238 }
1239 
1241 
1242  TFile *f;
1243  TString filename;
1244  if (!fwendcap){
1245  std::cout<< "File name = " << GetGeometryFileName().Data() << std::endl;
1246  filename = GetGeometryFileName();
1247  }else{
1248  std::cout<< "File name = " << fgeoName2 << std::endl;
1249  filename = fgeoName2;
1250  }
1251  f = new TFile(filename);
1252 
1253  TGeoVolume *Volume;
1254  TGeoCombiTrans *TransRotMatrix;
1255  if(filename.EndsWith("proto60.root")){
1256 
1257  f->GetObject("Proto60",Volume);
1258  TransRotMatrix = new TGeoCombiTrans(0.,0.,0.,new TGeoRotation());
1259  }else{
1260  Volume=(TGeoVolume *)f->Get("Emc3");
1261  TGeoRotation rotVolume;
1262  rotVolume.RotateY(180.);
1263  if(filename.Contains("proto")){
1264  rotVolume.RotateY(-6.51324892093148211e0);
1265  rotVolume.RotateX(5.15340666154607074e0);
1266  TransRotMatrix = new TGeoCombiTrans(-3.66092554252584108e+01, -2.84174235113817986e+01, 1.02907349180644744e+02,new TGeoRotation(rotVolume));
1267  //Moving Proto192 so that beam axis is on bottom right corner of row 49 crystal 21 (X4Y3-7)
1268 
1269  }else{
1270  TransRotMatrix = new TGeoCombiTrans(0., 0., 213.9,new TGeoRotation(rotVolume));//distance of the FwEndCap module to the target point was obtained to be around 2139 mm.}
1271  }
1272  }
1273  if(Volume == NULL){
1274  printf("Could not get geometry from file %s!.\nIs this the right file?\n",filename.Data());
1275  return;
1276  }
1277  TGeoVolume *Cave = gGeoManager->GetTopVolume();
1278  TGeoNode *n=Volume->GetNode(0);
1279 
1280  gGeoManager->AddVolume(Volume);
1281  TGeoVoxelFinder *voxels = Volume->GetVoxels();
1282  if (voxels) voxels->SetNeedRebuild();
1283  TGeoMatrix *M = n->GetMatrix();
1284  M->SetDefaultName();
1285  gGeoManager->GetListOfMatrices()->Remove(M);
1286  TGeoHMatrix *global = gGeoManager->GetHMatrix();
1287  gGeoManager->GetListOfMatrices()->Remove(global); //Remove the Identity matrix
1288 
1289 
1290  Cave->AddNode(Volume,0, TransRotMatrix);
1291 
1292  ExpandNode(Volume,Cave);
1293 
1294 }
1295 
1296 void PndEmc::ConstructRootGeomMod12() { // Added for barrel root file
1297  // volume from the root file
1298  TFile *f = new TFile(GetGeometryFileName());
1299  std::cout<< "File name = " << GetGeometryFileName().Data() << std::endl;
1300  TGeoVolume *Volume;
1301  Volume=(TGeoVolume *)f->Get("BarrelEMC");
1302  gGeoManager->AddVolume(Volume);
1303  TGeoVoxelFinder *voxels = Volume->GetVoxels();
1304  if (voxels) voxels->SetNeedRebuild();
1305 
1306  // cave
1307  TGeoVolume *Cave = gGeoManager->GetTopVolume();
1308  TGeoNode *n=Volume->GetNode(0);
1309  TGeoMatrix *M = n->GetMatrix();
1310  M->SetDefaultName();
1311  gGeoManager->GetListOfMatrices()->Remove(M);
1312  TGeoHMatrix *global = gGeoManager->GetHMatrix();
1313  gGeoManager->GetListOfMatrices()->Remove(global); //Remove the Identity matrix
1314 
1315  // place the volume to the cave
1316  TGeoCombiTrans* TransRotMatrix = new TGeoCombiTrans(0.,0.,0.,new TGeoRotation());
1317  Cave->AddNode(Volume,0, TransRotMatrix);
1318 
1319  ExpandNode(Volume,Cave);
1320 }
1321 
1323 
1324  TFile *fb;
1325  if (!bwendcap){
1326  std::cout<< "File name Bw = " << GetGeometryFileName().Data() << std::endl;
1327  fb=new TFile(GetGeometryFileName().Data());
1328  }else{
1329  std::cout<< "File name Bw1= " << fgeoName3 << std::endl;
1330  fb=new TFile(fgeoName3);
1331  }
1332 
1333  TGeoVolume *BwEmc=(TGeoVolume *)fb->Get("Emc4");
1334  TGeoVolume *Cave = gGeoManager->GetTopVolume();
1335  TGeoNode *n=BwEmc->GetNode(0);
1336 
1337 
1338 
1339  gGeoManager->AddVolume(BwEmc);
1340  TGeoVoxelFinder *voxels = BwEmc->GetVoxels();
1341  if (voxels) voxels->SetNeedRebuild();
1342  TGeoMatrix *M = n->GetMatrix();
1343  M->SetDefaultName();
1344  gGeoManager->GetListOfMatrices()->Remove(M);
1345  TGeoHMatrix *global = gGeoManager->GetHMatrix();
1346  gGeoManager->GetListOfMatrices()->Remove(global); //Remove the Identity matrix
1347 
1348  TGeoRotation rotBwEmc;
1349  rotBwEmc.RotateY(0.);
1350 
1351  //The first position of the BwEndCap crystals (center of the crystal!!) was -66 cm
1352  //Cave->AddNode(BwEmc,0, new TGeoCombiTrans(0., 0., -66.,new TGeoRotation(rotBwEmc)));
1353 
1354  // According to the last geometry integration of BwEndCap (22.04.09)
1355  // the position of the BwEndCap crystals (center of the crystal!!)
1356  // is -69.4 cm : -(3.+0.2+0.2+56.)-10 cm, which correspods to:
1357  // -(insulation+carbon fiber+safety distance+front face of crystals)-half size of crystal
1358  //--------------------------------------------------------------------
1359  // 2017 geometry version: the distance from the interaction point to the
1360  // outer surface of the bw EMC is fixed to 560 mm (could change in the future)
1361  static const Double_t distTargetBWEC = 56.0; // cm
1362  // The position of the Node is calculated from the z dimension of the outer volume
1363  TGeoNode *outerVol=BwEmc->FindNode("BWECouterVol_0");
1364  if(outerVol != NULL){
1365 
1366  TGeoCompositeShape* shape=(TGeoCompositeShape*)outerVol->GetVolume()->GetShape();
1367  Double_t origin[3]={0,0,0}, zdir[3]={0,0,1};
1368  Double_t halfLength = shape->DistFromInside(origin,zdir);
1369  Double_t bwEmcZpos = - (distTargetBWEC + halfLength);
1370  printf("PndEmc::ConstructRootGeomMod4: halfLength = %e, fullLength = %e, z-pos = %e, distTargetBWEC = %e\n",
1371  halfLength,2.*halfLength,bwEmcZpos,distTargetBWEC);
1372  Cave->AddNode(BwEmc,0, new TGeoCombiTrans(0., 0., bwEmcZpos,new TGeoRotation(rotBwEmc)));
1373 
1374  }else{ // for older geometry versions:
1375  Cave->AddNode(BwEmc,0, new TGeoCombiTrans(0., 0., -69.4,new TGeoRotation(rotBwEmc)));
1376  }
1377 
1378 
1379  ExpandNode(BwEmc,Cave);
1380 }
1381 
1383 
1384  TFile *fb;
1385  TString FileName = GetGeometryFileName().Data();
1386  if(FileName.EndsWith("5_fsc.root"))
1387  fb=new TFile(FileName);
1388  else if(fgeoName4.EndsWith("5_fsc.root")){
1389  fb=new TFile(fgeoName4);
1390  FileName = fgeoName4;
1391  }
1392 
1393  std::cout<< "File name Fsc= " << FileName << std::endl;
1394 
1395 
1396 
1397 
1398  TGeoVolume *Fsc=(TGeoVolume *)fb->Get("Emc5");
1399  TGeoVolume *Cave = gGeoManager->GetTopVolume();
1400  TGeoNode *n=Fsc->GetNode(0);
1401  if(fVerboseLevel>2) cout << "=====PndEmc::ConstructRootGeomMod5()====="<<endl;
1402  if(fVerboseLevel>2)Cave->Print();
1403  if(fVerboseLevel>2)Cave->PrintNodes();
1404 
1405  if(fVerboseLevel>2)Fsc->Print();
1406  if(fVerboseLevel>2)Fsc->PrintNodes();
1407 
1408  gGeoManager->AddVolume(Fsc);
1409  TGeoVoxelFinder *voxels = Fsc->GetVoxels();
1410  if (voxels) voxels->SetNeedRebuild();
1411  TGeoMatrix *M = n->GetMatrix();
1412  M->SetDefaultName();
1413  gGeoManager->GetListOfMatrices()->Remove(M);
1414  TGeoHMatrix *global = gGeoManager->GetHMatrix();
1415  gGeoManager->GetListOfMatrices()->Remove(global); //Remove the Identity matrix
1416 
1417  TGeoRotation rotFsc;
1418  Cave->AddNode(Fsc,0, new TGeoCombiTrans(0., 0., 818.775, new TGeoRotation(rotFsc)));
1419 
1420  ExpandNode(Fsc,Cave);
1421 
1422 }
1423 
1424 void PndEmc::ExpandNode(TGeoVolume *fVol, TGeoVolume *Cave){
1425 
1426  FairGeoLoader*geoLoad = FairGeoLoader::Instance();
1427  FairGeoInterface *geoFace = geoLoad->getGeoInterface();
1428  FairGeoMedia *Media = geoFace->getMedia();
1429  FairGeoBuilder *geobuild=geoLoad->getGeoBuilder();
1430 
1431  TObjArray *nodeList=fVol->GetNodes();
1432  Int_t nodes = nodeList->GetEntries();
1433 
1434  for (Int_t nod=0; nod < nodes; nod++){
1435 
1436  TGeoNode *fNode =(TGeoNode *)nodeList->At(nod);
1437  TGeoVolume *v= fNode->GetVolume();
1438 
1439  if(fNode->GetNdaughters()>0)
1440  ExpandNode(v,Cave);
1441 
1442  TGeoMedium* med1=v->GetMedium();
1443  //if(fVerboseLevel>2) std::cout<< "DEBUG NodeName = " << fNode->GetName() << std::endl;
1444  if (med1) {
1445  //if(fVerboseLevel>2) std::cout<< "DEBUG medium = " << med1->GetName() << std::endl;
1446  TGeoMaterial*mat1=v->GetMaterial();
1447  TGeoMaterial *newMat = gGeoManager->GetMaterial(mat1->GetName());
1448  if (newMat==0) {
1449  //std::cout<< "Material " << mat1->GetName() << " is not defined " << std::endl;
1450  FairGeoMedium *CbmMedium=Media->getMedium(mat1->GetName());
1451  if (!CbmMedium) {
1452  //std::cout << "Material is not defined in ASCII file nor in Root file" << std::endl;
1453  CbmMedium=new FairGeoMedium(mat1->GetName());
1454  Media->addMedium(CbmMedium);
1455  }
1456  //std::cout << "Create Medium " << mat1->GetName() << std::endl;
1457  Int_t nmed=geobuild->createMedium(CbmMedium);
1458  v->SetMedium(gGeoManager->GetMedium(nmed));
1459  gGeoManager->SetAllIndex();
1460  } else {
1461  //if(fVerboseLevel>2)
1462  //std::cout<< "DEBUG material was defined MaterialName= " << mat1->GetName() << std::endl;
1463  TGeoMedium *med2= gGeoManager->GetMedium(mat1->GetName());
1464  v->SetMedium(med2);
1465  }
1466  }
1467 
1468  TString name = v->GetName();
1469 
1470  if (!gGeoManager->FindVolumeFast(v->GetName())) {
1471  v->RegisterYourself();
1472  if(fVerboseLevel>2) std::cout<< "Volume -> " <<v->GetName()<< " --> Registered to gGeoManager" << endl;
1473  if (name.Contains("Crystal") || name.Contains("FscSci") || name.Contains("FscFiber")){
1474  AddSensitiveVolume(v);
1475  if(fVerboseLevel>2) std::cout<< "Volume " <<v->GetName()<< " is added to sensitives "<<endl;
1476  }
1477  }
1478  }
1479 }
1480 
1482  // Definition of materials
1483 
1484  FairGeoLoader*geoLoad = FairGeoLoader::Instance();
1485  FairGeoInterface *geoFace = geoLoad->getGeoInterface();
1486  FairGeoMedia *Media = geoFace->getMedia();
1487  FairGeoBuilder *geobuild=geoLoad->getGeoBuilder();
1488 
1489  FairGeoMedium *CbmMediumPb = Media->getMedium("lead");
1490  FairGeoMedium *CbmMediumPWO = Media->getMedium("PWO");
1491  FairGeoMedium *CbmMediumFsc = Media->getMedium("FscScint");
1492 
1493  //Int_t nmedPb= //[R.K.03/2017] unused variable
1494  geobuild->createMedium(CbmMediumPb);
1495  //Int_t nmedPWO= //[R.K.03/2017] unused variable
1496  geobuild->createMedium(CbmMediumPWO);
1497  //Int_t nmedFsc= //[R.K.03/2017] unused variable
1498  geobuild->createMedium(CbmMediumFsc);
1499 
1500  TGeoVolume *flayer1 = new TGeoVolumeAssembly("EmcLayer1");
1501  TGeoVolume *flayer2 = new TGeoVolumeAssembly("EmcLayer2");
1502  TGeoVolume *flayer2Hole = new TGeoVolumeAssembly("EmcLayer2Hole");
1503  TGeoVolume *flayer3 = new TGeoVolumeAssembly("Emc3");
1504  TGeoVolume *flayer4 = new TGeoVolumeAssembly("Emc4");
1505  TGeoVolume *flayer5 = new TGeoVolumeAssembly("Fsc");
1506  TGeoVolume *flayer6 = new TGeoVolumeAssembly("EmcTest");
1507 
1508  Bool_t bIsModuleOn[6] = {kFALSE, kFALSE, kFALSE, kFALSE, kFALSE,kFALSE};
1509  //Bool_t isFirst = kTRUE; //[R.K. 01/2017] unused variable?
1510 
1511  PndEmcReader read(GetGeometryFileName() );
1512  for(Int_t module=read.GetMinModules(); module<=read.GetMaxModules(); module++) {
1513  cout << "Emc module = " << module;
1514  if (module==5 && bIsFastFsc) {
1515  cout << "\t FAST" << endl << "******** " << endl;
1516  continue;
1517  }
1518  cout << endl << "******** " << endl;
1519 
1520  for(Int_t row=read.GetMinRows(module); row<=read.GetMaxRows(module); row++){
1521  for(Int_t crystal=read.GetMinCrystals(module,row); crystal<=read.GetMaxCrystals(module,row); crystal++) {
1522 
1523  Text_t buffer[30];
1524  sprintf(buffer,"emc0%dr%dc%d",module, row, crystal);
1525  DataG4 data = read.GetData(module,row,crystal);
1526 
1527  if (data.module==-1) continue; //if the pad is not present, do not create geometry
1528 
1529  if ((module<5) || (module==6) )
1530  { // Construction of target spectrometer geometry
1531  TGeoTrap *trap = new TGeoTrap(data.pDz/10., data.pTheta, data.pPhi,
1532  data.pDy1/10., data.pDx1/10., data.pDx2/10., data.pAlp1,
1533  data.pDy2/10., data.pDx3/10., data.pDx4/10., data.pAlp2);
1534  TGeoVolume *volume;
1535  volume = new TGeoVolume(buffer, trap, gGeoManager->GetMedium("PWO"));
1536 
1537  volume->SetLineColor(5);
1538  TGeoRotation rot;
1539  rot.RotateZ(data.tau);
1540  rot.RotateY(data.theta);
1541  rot.RotateZ(data.phi);
1542 
1543  if(module ==1) flayer1->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10.+3.7, new TGeoRotation (rot))); // shift of 37 mm with respect to interaction point
1544  if(module ==2) flayer2->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10.+3.7, new TGeoRotation (rot))); // shift of 37 mm with respect to interaction point
1545  if ((module ==2)&&
1546  !((crystal>=4 && crystal<=6) && (row<=3))
1547  )
1548  flayer2Hole->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10.+3.7, new TGeoRotation (rot))); // shift of 37 mm with respect to interaction point
1549  if(module ==3) flayer3->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10., new TGeoRotation (rot)));
1550  if(module ==4) flayer4->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10., new TGeoRotation (rot)));
1551  if(module ==6) flayer6->AddNode(volume,0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10., new TGeoRotation (rot)));
1552  bIsModuleOn[module-1] = kTRUE;
1553  AddSensitiveVolume(volume);
1554  }
1555 
1556  if (module==5 && !bIsFastFsc)
1557  { // Construction of forward spectrometer geometry
1558  // For the forward calorimeter the box geometry was written giving different meaning to the DataG4 members:
1559  // Dx1 Dy1 : half-widths in XY of the absorber
1560  // Dx2 Dy2 : half-widths of XY of the scintillator
1561  // pAlp1 pAlp2 : half-withs on Z of absorber/scintillator
1562  // pDz : full-width of space in between absorber and scintillator
1563  // pDx3 : number of layers of Pb/Scint
1564  TGeoVolume *volAbs, *volSci, *volCrystal;
1565  TGeoRotation rot;
1566  // TGeoVolume *volume = new TGeoVolumeAssembly(buffer);
1567  Double_t BoxZ=0.5*data.pDx3*(2.*data.pDz+2.*data.pAlp1+2.*data.pAlp2)/10.;
1568  TGeoBBox *volume = new TGeoBBox(data.pDx1/10., data.pDy1/10., BoxZ );
1569  volCrystal = new TGeoVolume(buffer, volume, gGeoManager->GetMedium("air"));
1570  volCrystal->SetLineColor(3);
1571  // cout << "\t Check " << buffer << "******** " << endl;
1572 
1573  TGeoBBox *absorber = new TGeoBBox(data.pDx1/10., data.pDy1/10., data.pAlp1/10.);
1574  TGeoBBox *scintillator = new TGeoBBox(data.pDx2/10., data.pDy2/10., data.pAlp2/10.);
1575  volAbs = new TGeoVolume("FscAbsorber", absorber, gGeoManager->GetMedium("lead"));
1576  volSci = new TGeoVolume("FscScintillator", scintillator, gGeoManager->GetMedium("FscScint"));
1577  volSci->SetLineColor(6);
1578 
1579  TGeoVolume *ffsclay = new TGeoVolumeAssembly("FscLayer");
1580  ffsclay->AddNode(volAbs, 0, new TGeoCombiTrans(0., 0., data.pAlp1/10., new TGeoRotation (rot)));
1581  ffsclay->AddNode(volSci, 0, new TGeoCombiTrans(0., 0., 2*data.pAlp1/10+data.pDz/10.+data.pAlp2/10., new TGeoRotation (rot)));
1582  AddSensitiveVolume(volSci);
1583 
1584  for (Int_t ll=0; ll<data.pDx3; ll++) {
1585  volCrystal->AddNode(ffsclay,ll, new TGeoCombiTrans(0., 0., ll*(2.*data.pDz+2.*data.pAlp1+2.*data.pAlp2)/10.-BoxZ, new TGeoRotation (rot)));
1586  //cout << " check " << volCrystal->GetName() << " has " << ffsclay->GetName() << " At " << ll*(2.*data.pDz+2.*data.pAlp1+2.*data.pAlp2)/10. << endl;
1587  }
1588  flayer5->AddNode(volCrystal, 0, new TGeoCombiTrans(data.posX/10., data.posY/10., data.posZ/10.+ BoxZ, new TGeoRotation (rot)));
1589  bIsModuleOn[module-1] = kTRUE;
1590  }
1591  }
1592  }
1593  }
1594 
1595  // Fast construction of the fsc geometry
1596  if (bIsFastFsc)
1597  {
1598  TGeoVolume *volAbs1, *volSci1, *volAbs2, *volSci2, *volAbs3, *volSci3;
1599  TGeoRotation rot;
1600  TGeoVolume *volume = new TGeoVolumeAssembly("FscBox");
1601  Int_t padX = 28,
1602  padY = 14,
1603  //holeX = 4, //[R.K. 01/2017] unused variable?
1604  holeY = 2,
1605  padX1 = 11,
1606  padX2 = 13;
1607  Float_t pDx1 = 11., pDx2 = 11., pDy1 = 11., pDy2 = 11., pAlp1 = 0.0275, pAlp2 = 0.15, pDz = 0.00375, posZ = 760.;
1608 
1609  TGeoBBox *absorber1 = new TGeoBBox(pDx1*padX/2. , pDy1*(padY-holeY)/4., pAlp1/2.);
1610  TGeoBBox *scintillator1 = new TGeoBBox(pDx2*padX/2. , pDy2*(padY-holeY)/4., pAlp2/2.);
1611  TGeoBBox *absorber2 = new TGeoBBox(pDx1*padX1/2., pDy1*holeY/4. , pAlp1/2.);
1612  TGeoBBox *scintillator2 = new TGeoBBox(pDx2*padX1/2., pDy2*holeY/4. , pAlp2/2.);
1613  TGeoBBox *absorber3 = new TGeoBBox(pDx1*padX2/2., pDy1*holeY/4. , pAlp1/2.);
1614  TGeoBBox *scintillator3 = new TGeoBBox(pDx2*padX2/2., pDy2*holeY/4. , pAlp2/2.);
1615 
1616  volAbs1 = new TGeoVolume("FscAbsorber1", absorber1, gGeoManager->GetMedium("lead"));
1617  volSci1 = new TGeoVolume("FscScintillator1", scintillator1, gGeoManager->GetMedium("FscScint"));
1618  volAbs2 = new TGeoVolume("FscAbsorber2", absorber2, gGeoManager->GetMedium("lead"));
1619  volSci2 = new TGeoVolume("FscScintillator2", scintillator2, gGeoManager->GetMedium("FscScint"));
1620  volAbs3 = new TGeoVolume("FscAbsorber3", absorber3, gGeoManager->GetMedium("lead"));
1621  volSci3 = new TGeoVolume("FscScintillator3", scintillator3, gGeoManager->GetMedium("FscScint"));
1622 
1623  //TGeoVolume *volSci1X = //[R.K.03/2017] unused variable
1624  gGeoManager->Division("FscScintillator1X" ,"FscScintillator1" , 1, padX , -pDx2*padX/2. , pDx2);
1625  TGeoVolume *volSci1XY = gGeoManager->Division("FscScintillator1XY","FscScintillator1X", 2, (padY-holeY)/2, -pDy2*(padY-holeY)/4., pDy2);
1626  //TGeoVolume *volSci2X = //[R.K.03/2017] unused variable
1627  gGeoManager->Division("FscScintillator2X" ,"FscScintillator2" , 1, padX1 , -pDx2*padX1/2. , pDx2);
1628  TGeoVolume *volSci2XY = gGeoManager->Division("FscScintillator2XY","FscScintillator2X", 2, holeY/2 , -pDy2*holeY/4. , pDy2);
1629  //TGeoVolume *volSci3X = //[R.K.03/2017] unused variable
1630  gGeoManager->Division("FscScintillator3X" ,"FscScintillator3" , 1, padX2 , -pDx2*padX2/2. , pDx2);
1631  TGeoVolume *volSci3XY = gGeoManager->Division("FscScintillator3XY","FscScintillator3X", 2, holeY/2 , -pDy2*holeY/4. , pDy2);
1632 
1633  TGeoVolume *ffsclay = new TGeoVolumeAssembly("FscLayer");
1634  ffsclay->AddNode(volAbs1, 0, new TGeoCombiTrans(0. , pDy1*(padY+holeY)/4., pAlp1/2. , new TGeoRotation (rot)));
1635  ffsclay->AddNode(volSci1, 0, new TGeoCombiTrans(0. , pDy2*(padY+holeY)/4., pAlp1+pDz+pAlp2/2., new TGeoRotation (rot)));
1636  ffsclay->AddNode(volAbs2, 0, new TGeoCombiTrans(pDx1*(-padX+padX1)/2., pDy1*holeY/4. , pAlp1/2. , new TGeoRotation (rot)));
1637  ffsclay->AddNode(volSci2, 0, new TGeoCombiTrans(pDx2*(-padX+padX1)/2., pDy2*holeY/4 , pAlp1+pDz+pAlp2/2., new TGeoRotation (rot)));
1638  ffsclay->AddNode(volAbs3, 0, new TGeoCombiTrans(pDx1*(padX-padX2)/2. , pDy1*holeY/4. , pAlp1/2. , new TGeoRotation (rot)));
1639  ffsclay->AddNode(volSci3, 0, new TGeoCombiTrans(pDx2*(padX-padX2)/2. , pDy2*holeY/4. , pAlp1+pDz+pAlp2/2., new TGeoRotation (rot)));
1640 
1641 
1642  AddSensitiveVolume(volSci1XY);
1643  AddSensitiveVolume(volSci2XY);
1644  AddSensitiveVolume(volSci3XY);
1645 
1646  for (Int_t ll=1; ll<=300; ll++) {
1647  volume->AddNode(ffsclay,ll, new TGeoCombiTrans(0., 0., ll*(2.*pDz+pAlp1+pAlp2), new TGeoRotation (rot)));
1648  }
1649 
1650  flayer5->AddNode(volume, 0, new TGeoCombiTrans(0., 0., posZ, new TGeoRotation (rot)));
1651  bIsModuleOn[4] = kTRUE;
1652  }
1653 
1654 
1655  TGeoVolume *flayer12 = new TGeoVolumeAssembly("Emc12");
1656  TGeoVolume *flayer12Hole = new TGeoVolumeAssembly("Emc12Hole");
1657  if (bIsModuleOn[0]) flayer12->AddNode(flayer1,0, new TGeoCombiTrans(0., 0., 0., new TGeoRotation(0)));
1658  if (bIsModuleOn[1]) flayer12->AddNode(flayer2,0, new TGeoCombiTrans(0., 0., 0., new TGeoRotation(0)));
1659  if (bIsModuleOn[0]) flayer12Hole->AddNode(flayer1,0, new TGeoCombiTrans(0., 0., 0., new TGeoRotation(0)));
1660  if (bIsModuleOn[1]) flayer12Hole->AddNode(flayer2Hole,0, new TGeoCombiTrans(0., 0., 0., new TGeoRotation(0)));
1661 
1662  TString vname = "cave";
1663  vname = vname.Strip();
1664  TGeoVolume* vcave = gGeoManager->FindVolumeFast(vname.Data());
1665  //if (bIsModuleOn[0] || bIsModuleOn[1]) vcave->AddNode(flayer12, 1);
1666  if (bIsModuleOn[2]) vcave->AddNode(flayer3, 1);
1667  if (bIsModuleOn[3]) vcave->AddNode(flayer4, 1);
1668  if (bIsModuleOn[4]) vcave->AddNode(flayer5, 1);
1669  if (bIsModuleOn[5]) vcave->AddNode(flayer6, 1);
1670 
1671  // 15 copies for barrel part of EMC (1st copy exists in emc_module12345.dat)
1672  if (bIsModuleOn[0] || bIsModuleOn[1])
1673  for (Int_t n=0;n<=15;n++){
1674  TGeoRotation rot1;
1675  rot1.RotateZ(22.5*n);
1676  if (n==0 || n==8) vcave->AddNode(flayer12Hole, n+1,new TGeoCombiTrans(0., 0., 0., new TGeoRotation (rot1)) );
1677  else vcave->AddNode(flayer12, n+1,new TGeoCombiTrans(0., 0., 0., new TGeoRotation (rot1)) );
1678  }
1679 
1680  // 3 copies for forward endcap of EMC (1st copy exists in emc_module12345.dat)
1681  if (bIsModuleOn[2])
1682  for (Int_t n=1;n<=3;n++){
1683  TGeoCombiTrans reflection;
1684  if (n==1) reflection.ReflectY(1);
1685  if (n==2) {reflection.ReflectY(1); reflection.ReflectX(1);}
1686  if (n==3) reflection.ReflectX(1);
1687  vcave->AddNode(flayer3, n+1, new TGeoCombiTrans(reflection));
1688  }
1689 
1690  // 3 copies for backward endcap of EMC (1st copy exists in emc_module12345.dat)
1691  if (bIsModuleOn[3])
1692  for (Int_t n=1;n<=3;n++){
1693  TGeoCombiTrans reflection;
1694  if (n==1) reflection.ReflectY(1);
1695  if (n==2) {reflection.ReflectY(1); reflection.ReflectX(1);}
1696  if (n==3) reflection.ReflectX(1);
1697  vcave->AddNode(flayer4, n+1, new TGeoCombiTrans(reflection));
1698  }
1699 
1700  // 3 copies for forward calorimeter (1st copy exists in emc_module12345.dat)
1701  if (bIsModuleOn[4])
1702  {
1703  if (bIsFastFsc)
1704  {
1705  TGeoCombiTrans reflection;
1706  reflection.ReflectY(1);
1707  vcave->AddNode(flayer5, 2, new TGeoCombiTrans(reflection));
1708  }
1709  else
1710  for (Int_t n=1;n<=3;n++){
1711  TGeoCombiTrans reflection;
1712  if (n==1) reflection.ReflectY(1);
1713  if (n==2) {reflection.ReflectY(1); reflection.ReflectX(1);}
1714  if (n==3) reflection.ReflectX(1);
1715  vcave->AddNode(flayer5, n+1, new TGeoCombiTrans(reflection));
1716  }
1717  }
1718 }
1719 
1720 // ----- Private method AddHit --------------------------------------------
1721 PndEmcPoint* PndEmc::AddHit(Int_t trackID, Int_t detID, Int_t evtID, TVector3 pos, TVector3 mom, Double_t time, Double_t length, Double_t eLoss, Short_t mod, Short_t row, Short_t crys, Short_t copy, Bool_t entering, Bool_t exiting) {
1722  TClonesArray& clref = *fEmcCollection;
1723  Int_t size = clref.GetEntriesFast();
1724  if (fVerboseLevel>1)
1725  cout << "-I- PndEmc: Adding Point at IN (" << pos.X() << ", " << pos.Y()
1726  << ", " << pos.Z() << ") cm, detector " << detID << ", evt " << evtID << ", track "
1727  << trackID <<", energy loss " << eLoss*1e06 << " keV, module " << mod << " row " << row << " crystal " << crys << " copy " << copy
1728  << ", entering " << entering << ", exiting " << exiting << endl;
1729 
1730  PndEmcPoint* myPoint = new(clref[size]) PndEmcPoint(trackID, detID, evtID, pos, mom, time, length, eLoss, mod, row, crys, copy, entering, exiting);
1731  myPoint->SetLink(FairLink(-1, FairRootManager::Instance()->GetEntryNr(), "MCTrack", trackID)); // 14.09.10 Stefano FIX
1732  return myPoint;
1733 }
1734 
1736  fTrackID = -999;
1737  fVolumeID = -999;
1738  fEventID = -999;
1739  fPos.SetXYZT(0., 0., 0., 0.);
1740  fMom.SetXYZT(0., 0., 0., 0.) ;
1741  fTime = -999;
1742  fLength = -999;
1743  fELoss = -999;
1744 }
1745 
1747  FairRun* fRun = FairRun::Instance();
1748 
1749  //check for GEANT3, else abort
1750  if (strcmp(fRun->GetName(),"TGeant3") == 0) {
1751 
1752  //get material ID for customs settings
1753 
1754  std::string sMedium[3] = {"FscScint", "FscFiber", "lead"};
1755  int matIdVMC;
1756  for(Int_t i = 0; i < 3; i++) {
1757  TGeoMedium *medium=gGeoManager->GetMedium(sMedium[i].c_str());
1758  if (medium==0) continue;
1759  matIdVMC = medium->GetId();
1760  double cut_el = 1.0E-4; // 100 KeV
1761  double cut_had = 1.0E-4; // 100 KeV
1762  //double tofmax = 1.E10; // (s) //[R.K. 01/2017] unused variable?
1763 
1764  // Set new properties, physics cuts etc. for the FSC
1765 // gMC->Gstpar(matIdVMC,"PAIR",1); /** pair production*/
1766 // gMC->Gstpar(matIdVMC,"COMP",1); /**Compton scattering*/
1767 // gMC->Gstpar(matIdVMC,"PHOT",1); /** photo electric effect */
1768 // gMC->Gstpar(matIdVMC,"PFIS",0); /**photofission*/
1769 // gMC->Gstpar(matIdVMC,"DRAY",1); /**delta-ray*/
1770 // gMC->Gstpar(matIdVMC,"ANNI",1); /**annihilation*/
1771 // gMC->Gstpar(matIdVMC,"BREM",1); /**bremsstrahlung*/
1772 // gMC->Gstpar(matIdVMC,"HADR",1); /**hadronic process*/
1773 // gMC->Gstpar(matIdVMC,"MUNU",1); /**muon nuclear interaction*/
1774 // gMC->Gstpar(matIdVMC,"DCAY",1); /**decay*/
1775 // gMC->Gstpar(matIdVMC,"LOSS",1); /**energy loss*/
1776 // gMC->Gstpar(matIdVMC,"MULS",1); /**multiple scattering*/
1777 // gMC->Gstpar(matIdVMC,"STRA",1); // collision sampling method to simulate energy loss in thin materials, particularly gases
1778 // gMC->Gstpar(matIdVMC,"RAYL",1); // Rayleigh scattering
1779 
1780  gMC->Gstpar(matIdVMC,"CUTGAM",cut_el);
1781  gMC->Gstpar(matIdVMC,"CUTELE",cut_el);
1782  gMC->Gstpar(matIdVMC,"CUTNEU",cut_had);
1783  gMC->Gstpar(matIdVMC,"CUTHAD",cut_had);
1784  gMC->Gstpar(matIdVMC,"CUTMUO",cut_el);
1785  gMC->Gstpar(matIdVMC,"BCUTE",cut_el);
1786  gMC->Gstpar(matIdVMC,"BCUTM",cut_el);
1787  gMC->Gstpar(matIdVMC,"DCUTE",cut_el);
1788  gMC->Gstpar(matIdVMC,"DCUTM",cut_el);
1789  gMC->Gstpar(matIdVMC,"PPCUTM",cut_el);
1790  }
1791  gMC->SetMaxNStep((int)1E6);
1792 
1793  std::cout<<"\n************************************************************\n"
1794  <<"PndEmc::SetSpecialPhysicsCuts():\n"
1795  <<" using special physics cuts ...\n";
1796  std::cout<<"************************************************************"
1797  <<std::endl;
1798 
1799  }
1800 }
1801 // ----
1802 
represents a mc hit in an emc crystal
Definition: PndEmcPoint.h:19
int row
Definition: anaLmdDigi.C:67
TVector3 pos
int GetMinRows(int module)
Bool_t bwendcap
Flag for the new FwEndCap geometry.
Definition: PndEmc.h:166
void ConstructRootGeomMod4()
Definition: PndEmc.cxx:1322
printf("RealTime=%f seconds, CpuTime=%f seconds\n", rtime, ctime)
TString fgeoName3
Definition: PndEmc.h:170
sim(Int_t nEvents=1, TString SimEngine="TGeant4", Float_t mom=6.231552)
double pDx3
Definition: PndEmcReader.h:22
double tau
Definition: PndEmcReader.h:20
FairGeoLoader * geoLoad
Int_t i
Definition: run_full.C:25
FairGeoMedia * Media
double pAlp2
Definition: PndEmcReader.h:22
void SetGeometryVersion(Int_t geometryVersion)
Definition: PndEmcGeoPar.h:26
Int_t fVolumeID
track index
Definition: PndEmc.h:153
Double_t val[nBoxes][nFEBox]
Definition: createCalib.C:11
void ConstructRootGeometry()
Definition: PndEmc.cxx:1240
virtual TClonesArray * GetCollection(Int_t iColl) const
Definition: PndEmc.cxx:910
virtual void ConstructGeometry()
Definition: PndEmc.cxx:1164
int col
Definition: anaLmdDigi.C:67
int GetMaxCrystals(int module, int row)
int GetMaxModules()
int n
double pDz
Definition: PndEmcReader.h:22
void SetMapperVersion(Int_t mapperVersion)
Definition: PndEmcGeoPar.h:21
Double_t par[3]
Double32_t fTime
momentum
Definition: PndEmc.h:157
double pDx2
Definition: PndEmcReader.h:22
Double_t mom
Definition: plot_dirc.C:14
void SetStorageOfData(Bool_t val)
Definition: PndEmc.cxx:941
TVector3 offset(2, 0, 0)
TString FileName
virtual void SetGeometryFileNameQuadruple(TString fname, TString fname2, TString fname3, TString fname4, TString geoVer="0")
Definition: PndEmc.cxx:1137
TGeoManager * gGeoManager
TString fgeoName4
Definition: PndEmc.h:171
__m128 v
Definition: P4_F32vec4.h:4
virtual void CopyClones(TClonesArray *cl1, TClonesArray *cl2, Int_t offset)
Definition: PndEmc.cxx:950
void AddPoint(DetectorId iDet)
Definition: PndStack.cxx:408
double phi
Definition: PndEmcReader.h:20
TGeoRotation * rot1
virtual void EndOfEvent()
Definition: PndEmc.cxx:895
Simulation of EMC.
Definition: PndEmc.h:26
virtual void Register()
Definition: PndEmc.cxx:902
int module
Definition: PndEmcReader.h:19
void ConstructRootGeomMod5()
Definition: PndEmc.cxx:1382
DataG4 GetData(int module, int row, int crystal)
FairRunAna * fRun
Definition: hit_dirc.C:58
Int_t fPosIndex
energy loss
Definition: PndEmc.h:160
TTree * T
Definition: anaLmdReco.C:32
Int_t a
Definition: anaLmdDigi.C:126
Bool_t bIsFastFsc
Hit collection.
Definition: PndEmc.h:163
FairGeoBuilder * geobuild
void ExpandNode(TGeoVolume *fVol, TGeoVolume *Cave)
Definition: PndEmc.cxx:1424
int nHits
Definition: RiemannTest.C:16
Cave SetGeometryFileName("pndcave.geo")
TGeoShape * shape
int GetMinCrystals(int module, int row)
double pPhi
Definition: PndEmcReader.h:22
int GetMaxRows(int module)
Double_t
TString medium
FairModule * Cave
Definition: sim_emc_apd.C:32
TGeoCombiTrans reflection
Int_t fTrackID
Definition: PndEmc.h:152
void ResetParameters()
Definition: PndEmc.cxx:1735
TLorentzVector fMom
position
Definition: PndEmc.h:156
virtual void SetGeometryVersion(const Int_t GeoNumber)
Definition: PndEmc.cxx:966
virtual void Initialize()
Definition: PndEmc.cxx:95
Double32_t fELoss
length
Definition: PndEmc.h:159
TString fgeoName2
Flag for the new BwEndCap geometry.
Definition: PndEmc.h:169
TFile * f
Definition: bump_analys.C:12
FairRuntimeDb * rtdb
Definition: hit_dirc.C:66
int getsign(const T &a)
Definition: PndEmc.cxx:53
virtual ~PndEmc()
Definition: PndEmc.cxx:83
double pDx4
Definition: PndEmcReader.h:22
TString name
double pDx1
Definition: PndEmcReader.h:22
int GetMinModules()
Int_t fEventID
volume id
Definition: PndEmc.h:154
double posX
Definition: PndEmcReader.h:21
virtual void SetTrackID(Int_t trackId)
Definition: PndEmcPoint.h:71
virtual void SetGeometryFileNameTriple(TString fname, TString fname2, TString fname3, TString geoVer="0")
Definition: PndEmc.cxx:1116
TGeoHMatrix * global
Bool_t fwendcap
Definition: PndEmc.h:165
double posZ
Definition: PndEmcReader.h:21
double pDy2
Definition: PndEmcReader.h:22
double pTheta
Definition: PndEmcReader.h:22
virtual void Print() const
Definition: PndEmc.cxx:918
TLorentzVector fPos
event id
Definition: PndEmc.h:155
virtual Bool_t ProcessHits(FairVolume *vol=0)
Definition: PndEmc.cxx:113
PndEmc()
Definition: PndEmc.cxx:59
virtual void BeginEvent()
Definition: PndEmc.cxx:106
ClassImp(PndAnaContFact)
TFile * fb
double pDy1
Definition: PndEmcReader.h:22
TGeoRotation rot
void ConstructRootGeomMod12()
Definition: PndEmc.cxx:1296
void ConstructASCIIGeometry()
Definition: PndEmc.cxx:1481
double posY
Definition: PndEmcReader.h:21
virtual void SetGeometryFileNameDouble(TString fname, TString fname2, Int_t fwbwchoice=0, TString geoVer="0")
Definition: PndEmc.cxx:1096
int sign(T val)
Definition: PndCADef.h:48
FairGeoMedium * CbmMediumPWO
FairGeoInterface * geoFace
virtual void Reset()
Definition: PndEmc.cxx:931
static int next[96]
Definition: ranlxd.cxx:374
TString nam
Definition: sim_hypGe.C:48
virtual void SetSpecialPhysicsCuts()
Definition: PndEmc.cxx:1746
Mvd Initialize()
double theta
Definition: PndEmcReader.h:20
double pAlp1
Definition: PndEmcReader.h:22
Int_t MapperVersion
Definition: PndEmc.h:173
PndEmcPoint * AddHit(Int_t trackID, Int_t detID, Int_t evtID, TVector3 pos, TVector3 mom, Double_t tof, Double_t length, Double_t eLoss, Short_t mod, Short_t row, Short_t crys, Short_t copy, Bool_t enterning, Bool_t exiting)
Definition: PndEmc.cxx:1721
TClonesArray * fEmcCollection
Definition: PndEmc.h:161
Double32_t fLength
time
Definition: PndEmc.h:158
const string filename
Bool_t fStoreData
Flag for fast fsc geometry.
Definition: PndEmc.h:164