SND@LHC Software
Loading...
Searching...
No Matches
splitcal Class Reference

#include <splitcal.h>

Inheritance diagram for splitcal:
Collaboration diagram for splitcal:

Public Member Functions

 splitcal (const char *Name, Bool_t Active)
 
 splitcal ()
 
virtual ~splitcal ()
 
virtual void Initialize ()
 
virtual Bool_t ProcessHits (FairVolume *v=0)
 
virtual void Register ()
 
virtual TClonesArray * GetCollection (Int_t iColl) const
 
virtual void Reset ()
 
void SetZStart (Double_t ZStart)
 
void SetEmpty (Double_t Empty, Double_t BigGap, Double_t ActiveECAL_gas_gap, Double_t first_precision_layer, Double_t second_precision_layer, Double_t third_precision_layer, Double_t num_precision_layers)
 
void SetThickness (Double_t ActiveECALThickness, Double_t ActiveHCALThickness, Double_t FilterECALThickness, Double_t FilterECALThickness_first, Double_t FilterHCALThickness, Double_t ActiveECAL_gas_Thickness)
 
void SetMaterial (Double_t ActiveECALMaterial, Double_t ActiveHCALMaterial, Double_t FilterECALMaterial, Double_t FilterHCALMaterial)
 
void SetNSamplings (Int_t nECALSamplings, Int_t nHCALSamplings, Double_t ActiveHCAL)
 
void SetNModules (Int_t nModulesInX, Int_t nModulesInY)
 
void SetNStrips (Int_t nStrips)
 
void SetStripSize (Double_t stripHalfWidth, Double_t stripHalfLength)
 
void SetXMax (Double_t xMax)
 
void SetYMax (Double_t yMax)
 
void ConstructGeometry ()
 
splitcalPointAddHit (Int_t trackID, Int_t detID, TVector3 pos, TVector3 mom, Double_t time, Double_t length, Double_t eLoss, Int_t pdgcode)
 
virtual void CopyClones (TClonesArray *cl1, TClonesArray *cl2, Int_t offset)
 
virtual void SetSpecialPhysicsCuts ()
 
virtual void EndOfEvent ()
 
virtual void FinishPrimary ()
 
virtual void FinishRun ()
 
virtual void BeginPrimary ()
 
virtual void PostTrack ()
 
virtual void PreTrack ()
 
virtual void BeginEvent ()
 

Private Member Functions

 splitcal (const splitcal &)
 
splitcaloperator= (const splitcal &)
 
Int_t InitMedium (const char *name)
 

Private Attributes

Int_t fTrackID
 
Int_t fVolumeID
 track index
 
TLorentzVector fPos
 volume id
 
TLorentzVector fMom
 position at entrance
 
Double_t fTime
 momentum at entrance
 
Double_t fLength
 time
 
Double_t fELoss
 length
 
Double_t fActiveECALThickness
 energy loss
 
Double_t fActiveHCALThickness
 
Double_t fFilterECALThickness
 
Double_t xfFilterECALThickness
 
Double_t fFilterECALThickness_first
 
Double_t fFilterHCALThickness
 
Double_t fActiveECALMaterial
 
Double_t fActiveHCALMaterial
 
Double_t fFilterECALMaterial
 
Double_t fFilterHCALMaterial
 
Double_t fActiveECAL_gas_gap
 
Double_t fActiveECAL_gas_Thickness
 
Int_t fnECALSamplings
 
Int_t fnHCALSamplings
 
Double_t fActiveHCAL
 
Double_t fZStart
 
Double_t fEmpty
 
Double_t fBigGap
 
Double_t fXMax
 
Double_t fYMax
 
Double_t ffirst_precision_layer
 
Double_t fsecond_precision_layer
 
Double_t fthird_precision_layer
 
Double_t fnum_precision_layers
 
Int_t fNModulesInX
 
Int_t fNModulesInY
 
Int_t fNStripsPerModule
 
Double_t fStripHalfWidth
 
Double_t fStripHalfLength
 
TClonesArray * fsplitcalPointCollection
 

Detailed Description

Definition at line 13 of file splitcal.h.

Constructor & Destructor Documentation

◆ splitcal() [1/3]

splitcal::splitcal ( const char *  Name,
Bool_t  Active 
)

Name : Detector Name Active: kTRUE for active detectors (ProcessHits() will be called) kFALSE for inactive detectors

Definition at line 52 of file splitcal.cxx.

53 : FairDetector(name, active, kSplitCal),
54 fTrackID(-1),
55 fVolumeID(-1),
56 fPos(),
57 fMom(),
58 fTime(-1.),
59 fLength(-1.),
60 fELoss(-1),
61 fsplitcalPointCollection(new TClonesArray("splitcalPoint"))
62{
63}
@ kSplitCal
Double_t fLength
time
Definition splitcal.h:107
TLorentzVector fMom
position at entrance
Definition splitcal.h:105
TLorentzVector fPos
volume id
Definition splitcal.h:104
Double_t fELoss
length
Definition splitcal.h:108
Int_t fTrackID
Definition splitcal.h:102
Double_t fTime
momentum at entrance
Definition splitcal.h:106
TClonesArray * fsplitcalPointCollection
Definition splitcal.h:125
Int_t fVolumeID
track index
Definition splitcal.h:103

◆ splitcal() [2/3]

splitcal::splitcal ( )

default constructor

Definition at line 39 of file splitcal.cxx.

40 : FairDetector("splitcal", kTRUE, kSplitCal),
41 fTrackID(-1),
42 fVolumeID(-1),
43 fPos(),
44 fMom(),
45 fTime(-1.),
46 fLength(-1.),
47 fELoss(-1),
48 fsplitcalPointCollection(new TClonesArray("splitcalPoint"))
49{
50}

◆ ~splitcal()

splitcal::~splitcal ( )
virtual

destructor

Definition at line 65 of file splitcal.cxx.

66{
70 }
71}

◆ splitcal() [3/3]

splitcal::splitcal ( const splitcal )
private

Member Function Documentation

◆ AddHit()

splitcalPoint * splitcal::AddHit ( Int_t  trackID,
Int_t  detID,
TVector3  pos,
TVector3  mom,
Double_t  time,
Double_t  length,
Double_t  eLoss,
Int_t  pdgcode 
)

This method is an example of how to add your own point of type splitcalPoint to the clones array

Definition at line 456 of file splitcal.cxx.

460{
461 TClonesArray& clref = *fsplitcalPointCollection;
462 Int_t size = clref.GetEntriesFast();
463 // cout << "splitcal hit called"<< pos.z()<<endl;
464 return new(clref[size]) splitcalPoint(trackID, detID, pos, mom,
465 time, length, eLoss, pdgCode);
466}

◆ BeginEvent()

virtual void splitcal::BeginEvent ( )
inlinevirtual

Definition at line 94 of file splitcal.h.

94{;}

◆ BeginPrimary()

virtual void splitcal::BeginPrimary ( )
inlinevirtual

Definition at line 91 of file splitcal.h.

91{;}

◆ ConstructGeometry()

void splitcal::ConstructGeometry ( )

Create the detector geometry

If you are using the standard ASCII input for the geometry just copy this and use it for your detector, otherwise you can implement here you own way of constructing the geometry.

Definition at line 254 of file splitcal.cxx.

255{
260 TGeoVolume *top=gGeoManager->GetTopVolume();
261 TGeoVolume *tSplitCal = new TGeoVolumeAssembly("SplitCalDetector");
262
263 InitMedium("iron");
264 InitMedium("lead");
265 InitMedium("Scintillator");
266 InitMedium("argon");
267 InitMedium("GEMmixture");
268
269 TGeoMedium *A2 =gGeoManager->GetMedium("iron");
270 TGeoMedium *A3 =gGeoManager->GetMedium("lead");
271 TGeoMedium *A4 =gGeoManager->GetMedium("GEMmixture");
272 TGeoMedium *A1 =gGeoManager->GetMedium("Scintillator");
273
274
275
276 Double_t zStartSplitCal = fZStart;
277
278 TGeoVolume *newECALfilter_first; //first layer can have different thickeness
279 TGeoVolume *newECALfilter;
280 TGeoVolume *newECALdet_gas;
281 TGeoVolume *stripGivingX;
282 TGeoVolume *stripGivingY;
283
284 TGeoVolume *newHCALfilter[100];
285 TGeoVolume *newHCALdet[100];
286 const char* char_labelHCALfilter[100];
287 TString labelHCALfilter;
288 const char* char_labelHCALdet[100];
289 TString labelHCALdet;
290
291
292 Double_t z_splitcal=0;
293 Int_t i_nlayECAL_gas;
294
295 // logical volume for the absorbing layers
296 // first absorbing layer can have different thinkens from the others
297 newECALfilter_first = gGeoManager->MakeBox("ECALfilter_first", A3, fXMax, fYMax, fFilterECALThickness_first/2);
298 newECALfilter_first->SetLineColor(kGray);
299 newECALfilter = gGeoManager->MakeBox("ECALfilter", A3, fXMax, fYMax, fFilterECALThickness/2);
300 newECALfilter->SetLineColor(kGray);
301
302 stripGivingX = gGeoManager->MakeBox("stripGivingX", A1, fStripHalfWidth, fStripHalfLength, fActiveECALThickness/2);
303 stripGivingX->SetVisibility(kTRUE);
304 AddSensitiveVolume(stripGivingX);
305 stripGivingX->SetLineColor(kGreen);
306
307 stripGivingY = gGeoManager->MakeBox("stripGivingY", A1, fStripHalfLength, fStripHalfWidth, fActiveECALThickness/2);
308 stripGivingY->SetVisibility(kTRUE);
309 AddSensitiveVolume(stripGivingY);
310 stripGivingY->SetLineColor(kGreen);
311
312 // logical volume for the high precision sensitive layers
313 newECALdet_gas = gGeoManager->MakeBox("ECALdet_gas", A4, fXMax, fYMax, fActiveECAL_gas_Thickness/2);
314 AddSensitiveVolume(newECALdet_gas);
315 newECALdet_gas->SetLineColor(kRed);
316
317
318 // now position layer volumes in tSplitCal mother volume
319 for (Int_t i_nlayECAL=0; i_nlayECAL<fnECALSamplings; i_nlayECAL++){
320
321 // position absorber layers
322 // thinkness of first layer can be different from others
323 if (i_nlayECAL==0){
324 z_splitcal+=fFilterECALThickness_first/2;
325 tSplitCal->AddNode(newECALfilter_first, i_nlayECAL*1e5, new TGeoTranslation(0, 0, z_splitcal));
326 z_splitcal+=fFilterECALThickness_first/2;
327 } else {
328 z_splitcal+=fFilterECALThickness/2;
329 tSplitCal->AddNode(newECALfilter, i_nlayECAL*1e5, new TGeoTranslation(0, 0, z_splitcal));
330 z_splitcal+=fFilterECALThickness/2;
331 }
332 // std::cout<< "--- i_nlayECAL*1e5 = "<< i_nlayECAL*1e5 << std::endl;
333
334 if(i_nlayECAL==0) z_splitcal+=fEmpty; // space after first layer? set to 0 in the config file? for whar is it for?
335 if(i_nlayECAL==7) z_splitcal+=fBigGap;
336
337 // position high precision sensitive layers
338 i_nlayECAL_gas=-1000;
339 if(i_nlayECAL==ffirst_precision_layer or i_nlayECAL==fsecond_precision_layer or i_nlayECAL==fthird_precision_layer){
340 z_splitcal+=fActiveECAL_gas_Thickness/2;
341 if(i_nlayECAL==ffirst_precision_layer) i_nlayECAL_gas=0;
342 else if(i_nlayECAL==fsecond_precision_layer ){
343 if(fnum_precision_layers==2) i_nlayECAL_gas=3;
344 else i_nlayECAL_gas=1;
345 }
346 else if(i_nlayECAL==fthird_precision_layer){
347 if(fnum_precision_layers==2) i_nlayECAL_gas=4;
348 else i_nlayECAL_gas=2;
349 }
350
351 tSplitCal->AddNode(newECALdet_gas, 1e8+(i_nlayECAL+1)*1e5 , new TGeoTranslation(0, 0, z_splitcal));
352 z_splitcal+=fActiveECAL_gas_Thickness/2;
354 z_splitcal+=fActiveECAL_gas_gap;
355 z_splitcal+=fActiveECAL_gas_Thickness/2;
356 if(i_nlayECAL==ffirst_precision_layer) i_nlayECAL_gas=1;
357 if(i_nlayECAL==fsecond_precision_layer) i_nlayECAL_gas=3;
358 tSplitCal->AddNode(newECALdet_gas, 1e8+(i_nlayECAL+1)*1e5 , new TGeoTranslation(0, 0, z_splitcal));
359 z_splitcal+=fActiveECAL_gas_Thickness/2;
360 }
361 // std::cout<< "--- index high precision layer = "<< 1e8+(i_nlayECAL+1)*1e5 << std::endl;
362 }
363 else{
364 // position sensitive layers
365 z_splitcal+=fActiveECALThickness/2;
366 if (i_nlayECAL%2==0) {
367 //strips giving x information
368 for (int mx=0; mx<fNModulesInX; mx++){
369 for(int my=0; my<fNModulesInY; my++){
370 for(int j=0;j<fNStripsPerModule;j++){
371 int index = (i_nlayECAL+1)*1e5 + (mx+1)*1e4 + (my+1)*1e3 + j+1 ;
372 double xCoordinate = -fXMax + (fNStripsPerModule*mx + j+0.5) * fStripHalfWidth * 2; // the times 2 is to get the total width from the half-width
373 double yCoordinate = -fYMax + (my+0.5) * fStripHalfLength * 2; // the times 2 is to get the total length from the half-length
374 tSplitCal->AddNode(stripGivingX, index, new TGeoTranslation(xCoordinate, yCoordinate, z_splitcal));
375 // std::cout<< "--- index = "<< index << std::endl;
376 }//end loop on strips
377 }//end loop on modules in y
378 }//end loop on modules in x
379 }//end layer stripped in X
380 else {
381 //strips giving y information
382 for (int mx=0; mx<fNModulesInX; mx++){
383 for(int my=0; my<fNModulesInY; my++){
384 for(int j=0;j<fNStripsPerModule;j++){
385 int index = (i_nlayECAL+1)*1e5 + (mx+1)*1e4 + (my+1)*1e3 + j+1 ;
386 double xCoordinate = -fXMax + (mx+0.5) * fStripHalfLength * 2; // the times 2 is to get the total length from the half-length
387 double yCoordinate = -fYMax + (fNStripsPerModule*my+ j+0.5) * fStripHalfWidth * 2; // the times 2 is to get the total width from the half-width
388 tSplitCal->AddNode(stripGivingY, index, new TGeoTranslation(xCoordinate, yCoordinate, z_splitcal));
389 // std::cout<< "--- index = "<< index << std::endl;
390 }//end loop on strips
391 }//end loop on modules in y
392 }//end loop on modules in x
393 }//end layer stripped in Y
394
395 z_splitcal+=fActiveECALThickness/2;
396 }// end loop on sensitive scintillator layers
397
398 } //end loop in ecal layers
399
400 for (Int_t i_nlayHCAL=0; i_nlayHCAL<1; i_nlayHCAL++){
401 labelHCALfilter="HCALfilter_";
402 labelHCALfilter+=i_nlayHCAL;
403 char_labelHCALfilter[i_nlayHCAL]=labelHCALfilter;
404 labelHCALdet="HCAL_det";
405 labelHCALdet+=i_nlayHCAL;
406 char_labelHCALdet[i_nlayHCAL]=labelHCALdet;
407 newHCALfilter[i_nlayHCAL] = gGeoManager->MakeBox(char_labelHCALfilter[i_nlayHCAL], A2, fXMax, fYMax, fFilterHCALThickness/2);
408 if(fActiveHCAL){
409 newHCALdet[i_nlayHCAL] = gGeoManager->MakeBox(char_labelHCALdet[i_nlayHCAL], A4, fXMax, fYMax, fActiveHCALThickness/2);
410
411 AddSensitiveVolume(newHCALdet[i_nlayHCAL]);
412
413 newHCALdet[i_nlayHCAL]->SetLineColor(kRed);
414 }
415 newHCALfilter[i_nlayHCAL]->SetLineColor(kBlue);
416 }
417 for (Int_t i_nlayHCAL=0; i_nlayHCAL<fnHCALSamplings; i_nlayHCAL++){
418 z_splitcal+=fFilterHCALThickness/2;
419 tSplitCal->AddNode(newHCALfilter[i_nlayHCAL], 1, new TGeoTranslation(0, 0, z_splitcal));
420 z_splitcal+=fFilterHCALThickness/2;
421 // z_splitcal+=fEmpty;
422 z_splitcal+=fActiveHCALThickness/2;
423 if(fActiveHCAL) tSplitCal->AddNode(newHCALdet[i_nlayHCAL], 1, new TGeoTranslation(0, 0, z_splitcal));
424 z_splitcal+=fActiveHCALThickness/2;
425 // z_splitcal+=fEmpty;
426
427
428 }
429
430
431
432 //finish assembly and position
433 TGeoShapeAssembly* asmb = dynamic_cast<TGeoShapeAssembly*>(tSplitCal->GetShape());
434 Double_t totLength = asmb->GetDZ();
435 top->AddNode(tSplitCal, 1, new TGeoTranslation(0, 0,zStartSplitCal+totLength));
436
437
438
439
440 // //gROOT->SetBatch(true);
441
442 // TCanvas* c1 = new TCanvas("splitcalCanvas", "", 800, 800);
443 // c1->cd();
444
445 // TView3D* tview = (TView3D*) TView::CreateView();
446 // tview->SetRange(-fXMax*1.2, -fYMax*1.2, 2500, fXMax*1.2, fYMax*1.2, 3800);
447 // tview->RotateView(0, 90, c1);
448
449 // tSplitCal->Draw("ogle");
450 // c1->SaveAs(TString("splitcal.eps"));
451 // c1->SaveAs(TString("splitcal.pdf"));
452
453
454}
Double_t fFilterHCALThickness
Definition splitcal.h:109
Double_t fnum_precision_layers
Definition splitcal.h:118
Double_t fActiveHCALThickness
Definition splitcal.h:109
Int_t fNModulesInY
Definition splitcal.h:119
Double_t fFilterECALThickness_first
Definition splitcal.h:109
Double_t fActiveECALThickness
energy loss
Definition splitcal.h:109
Double_t fthird_precision_layer
Definition splitcal.h:118
Double_t fActiveHCAL
Definition splitcal.h:113
Int_t fNStripsPerModule
Definition splitcal.h:120
Int_t fNModulesInX
Definition splitcal.h:119
Double_t fActiveECAL_gas_gap
Definition splitcal.h:111
Double_t fBigGap
Definition splitcal.h:115
Double_t fActiveECAL_gas_Thickness
Definition splitcal.h:111
Double_t fStripHalfLength
Definition splitcal.h:121
Double_t fEmpty
Definition splitcal.h:115
Double_t fZStart
Definition splitcal.h:114
Double_t fYMax
Definition splitcal.h:117
Int_t fnHCALSamplings
Definition splitcal.h:112
Double_t fXMax
Definition splitcal.h:116
Int_t InitMedium(const char *name)
Definition splitcal.cxx:80
Double_t fsecond_precision_layer
Definition splitcal.h:118
Double_t ffirst_precision_layer
Definition splitcal.h:118
Double_t fStripHalfWidth
Definition splitcal.h:121
Double_t fFilterECALThickness
Definition splitcal.h:109
Int_t fnECALSamplings
Definition splitcal.h:112

◆ CopyClones()

virtual void splitcal::CopyClones ( TClonesArray *  cl1,
TClonesArray *  cl2,
Int_t  offset 
)
inlinevirtual

The following methods can be implemented if you need to make any optional action in your detector during the transport.

Definition at line 85 of file splitcal.h.

86 {;}

◆ EndOfEvent()

void splitcal::EndOfEvent ( )
virtual

Definition at line 149 of file splitcal.cxx.

150{
151
153
154}

◆ FinishPrimary()

virtual void splitcal::FinishPrimary ( )
inlinevirtual

Definition at line 89 of file splitcal.h.

89{;}

◆ FinishRun()

virtual void splitcal::FinishRun ( )
inlinevirtual

Definition at line 90 of file splitcal.h.

90{;}

◆ GetCollection()

TClonesArray * splitcal::GetCollection ( Int_t  iColl) const
virtual

Gets the produced collections

Definition at line 173 of file splitcal.cxx.

174{
175 if (iColl == 0) { return fsplitcalPointCollection; }
176 else { return NULL; }
177}

◆ Initialize()

void splitcal::Initialize ( )
virtual

Initialization of the detector is done here

Definition at line 73 of file splitcal.cxx.

74{
75 FairDetector::Initialize();
76// FairRuntimeDb* rtdb= FairRun::Instance()->GetRuntimeDb();
77// splitcalGeoPar* par=(splitcalGeoPar*)(rtdb->getContainer("splitcalGeoPar"));
78}

◆ InitMedium()

Int_t splitcal::InitMedium ( const char *  name)
private

Definition at line 80 of file splitcal.cxx.

81{
82 static FairGeoLoader *geoLoad=FairGeoLoader::Instance();
83 static FairGeoInterface *geoFace=geoLoad->getGeoInterface();
84 static FairGeoMedia *media=geoFace->getMedia();
85 static FairGeoBuilder *geoBuild=geoLoad->getGeoBuilder();
86
87 FairGeoMedium *ShipMedium=media->getMedium(name);
88
89 if (!ShipMedium)
90 {
91 Fatal("InitMedium","Material %s not defined in media file.", name);
92 return -1111;
93 }
94 TGeoMedium* medium=gGeoManager->GetMedium(name);
95 if (medium!=NULL)
96 return ShipMedium->getMediumIndex();
97
98 return geoBuild->createMedium(ShipMedium);
99}

◆ operator=()

splitcal & splitcal::operator= ( const splitcal )
private

◆ PostTrack()

virtual void splitcal::PostTrack ( )
inlinevirtual

Definition at line 92 of file splitcal.h.

92{;}

◆ PreTrack()

virtual void splitcal::PreTrack ( )
inlinevirtual

Definition at line 93 of file splitcal.h.

93{;}

◆ ProcessHits()

Bool_t splitcal::ProcessHits ( FairVolume *  v = 0)
virtual

this method is called for each step during simulation (see FairMCApplication::Stepping())

This method is called from the MC stepping

fVolumeID = gGeoManager->FindVolumeFast(vol->GetName())->GetNumber();

Definition at line 100 of file splitcal.cxx.

101{
103 //Set parameters at entrance of volume. Reset ELoss.
104 if ( gMC->IsTrackEntering() ) {
105 fELoss = 0.;
106 fTime = gMC->TrackTime() * 1.0e09;
107 fLength = gMC->TrackLength();
108 gMC->TrackPosition(fPos);
109 gMC->TrackMomentum(fMom);
110 }
111
112 // Sum energy loss for all steps in the active volume
113 fELoss += gMC->Edep();
114
115 // Create splitcalPoint at exit of active volume
116 if ( gMC->IsTrackExiting() ||
117 gMC->IsTrackStop() ||
118 gMC->IsTrackDisappeared() ) {
119 fTrackID = gMC->GetStack()->GetCurrentTrackNumber();
120 //fVolumeID = vol->getMCid();
121 //cout << "splitcal proc "<< fVolumeID<<" "<<vol->GetName()<<" "<<vol->getVolumeId() <<endl;
122 // cout << " "<< gGeoManager->FindVolumeFast(vol->GetName())->GetNumber()<< " " << gMC->CurrentVolID() << endl;
124 //VolumeID = vol->getMCid();
125 Int_t detID=0;
126 gMC->CurrentVolID(detID);
127
128 //if (fVolumeID == detID) {
129 // return kTRUE; }
130 fVolumeID = detID;
131 // cout << " "<<fVolumeID << endl;
132 if (fELoss == 0. ) { return kFALSE; }
133 TParticle* p=gMC->GetStack()->GetCurrentTrack();
134 Int_t pdgCode = p->GetPdgCode();
135 // if(fVolumeID<405 && fTime<70){
136 AddHit(fTrackID, fVolumeID, TVector3(fPos.X(), fPos.Y(), fPos.Z()),
137 TVector3(fMom.Px(), fMom.Py(), fMom.Pz()), fTime, fLength,
138 fELoss,pdgCode);
139
140 // Increment number of splitcal det points in TParticle
141 ShipStack* stack = (ShipStack*) gMC->GetStack();
142 stack->AddPoint(kSplitCal);
143
144 }
145
146 return kTRUE;
147}
splitcalPoint * AddHit(Int_t trackID, Int_t detID, TVector3 pos, TVector3 mom, Double_t time, Double_t length, Double_t eLoss, Int_t pdgcode)
Definition splitcal.cxx:456

◆ Register()

void splitcal::Register ( )
virtual

Registers the produced collections in FAIRRootManager.

This will create a branch in the output tree called splitcalPoint, setting the last parameter to kFALSE means: this collection will not be written to the file, it will exist only during the simulation.

Definition at line 158 of file splitcal.cxx.

159{
160
167 FairRootManager::Instance()->Register("splitcalPoint", "splitcal",
169
170}

◆ Reset()

void splitcal::Reset ( )
virtual

has to be called after each event to reset the containers

Definition at line 179 of file splitcal.cxx.

180{
182}

◆ SetEmpty()

void splitcal::SetEmpty ( Double_t  Empty,
Double_t  BigGap,
Double_t  ActiveECAL_gas_gap,
Double_t  first_precision_layer,
Double_t  second_precision_layer,
Double_t  third_precision_layer,
Double_t  num_precision_layers 
)

◆ SetMaterial()

void splitcal::SetMaterial ( Double_t  ActiveECALMaterial,
Double_t  ActiveHCALMaterial,
Double_t  FilterECALMaterial,
Double_t  FilterHCALMaterial 
)

◆ SetNModules()

void splitcal::SetNModules ( Int_t  nModulesInX,
Int_t  nModulesInY 
)

Definition at line 229 of file splitcal.cxx.

230{
231 fNModulesInX=nModulesInX;
232 fNModulesInY=nModulesInY;
233}

◆ SetNSamplings()

void splitcal::SetNSamplings ( Int_t  nECALSamplings,
Int_t  nHCALSamplings,
Double_t  ActiveHCAL 
)

◆ SetNStrips()

void splitcal::SetNStrips ( Int_t  nStrips)

Definition at line 235 of file splitcal.cxx.

236{
237 fNStripsPerModule=nStrips;
238}

◆ SetSpecialPhysicsCuts()

virtual void splitcal::SetSpecialPhysicsCuts ( )
inlinevirtual

Definition at line 87 of file splitcal.h.

87{;}

◆ SetStripSize()

void splitcal::SetStripSize ( Double_t  stripHalfWidth,
Double_t  stripHalfLength 
)

Definition at line 240 of file splitcal.cxx.

241{
242 fStripHalfWidth=stripHalfWidth;
243 fStripHalfLength=stripHalfLength;
244}

◆ SetThickness()

void splitcal::SetThickness ( Double_t  ActiveECALThickness,
Double_t  ActiveHCALThickness,
Double_t  FilterECALThickness,
Double_t  FilterECALThickness_first,
Double_t  FilterHCALThickness,
Double_t  ActiveECAL_gas_Thickness 
)

◆ SetXMax()

void splitcal::SetXMax ( Double_t  xMax)

Definition at line 246 of file splitcal.cxx.

247{
248 fXMax=xMax;
249}

◆ SetYMax()

void splitcal::SetYMax ( Double_t  yMax)

Definition at line 250 of file splitcal.cxx.

251{
252 fYMax=yMax;
253}

◆ SetZStart()

void splitcal::SetZStart ( Double_t  ZStart)

Definition at line 184 of file splitcal.cxx.

Member Data Documentation

◆ fActiveECAL_gas_gap

Double_t splitcal::fActiveECAL_gas_gap
private

Definition at line 111 of file splitcal.h.

◆ fActiveECAL_gas_Thickness

Double_t splitcal::fActiveECAL_gas_Thickness
private

Definition at line 111 of file splitcal.h.

◆ fActiveECALMaterial

Double_t splitcal::fActiveECALMaterial
private

Definition at line 110 of file splitcal.h.

◆ fActiveECALThickness

Double_t splitcal::fActiveECALThickness
private

energy loss

Definition at line 109 of file splitcal.h.

◆ fActiveHCAL

Double_t splitcal::fActiveHCAL
private

Definition at line 113 of file splitcal.h.

◆ fActiveHCALMaterial

Double_t splitcal::fActiveHCALMaterial
private

Definition at line 110 of file splitcal.h.

◆ fActiveHCALThickness

Double_t splitcal::fActiveHCALThickness
private

Definition at line 109 of file splitcal.h.

◆ fBigGap

Double_t splitcal::fBigGap
private

Definition at line 115 of file splitcal.h.

◆ fELoss

Double_t splitcal::fELoss
private

length

Definition at line 108 of file splitcal.h.

◆ fEmpty

Double_t splitcal::fEmpty
private

Definition at line 115 of file splitcal.h.

◆ fFilterECALMaterial

Double_t splitcal::fFilterECALMaterial
private

Definition at line 110 of file splitcal.h.

◆ fFilterECALThickness

Double_t splitcal::fFilterECALThickness
private

Definition at line 109 of file splitcal.h.

◆ fFilterECALThickness_first

Double_t splitcal::fFilterECALThickness_first
private

Definition at line 109 of file splitcal.h.

◆ fFilterHCALMaterial

Double_t splitcal::fFilterHCALMaterial
private

Definition at line 110 of file splitcal.h.

◆ fFilterHCALThickness

Double_t splitcal::fFilterHCALThickness
private

Definition at line 109 of file splitcal.h.

◆ ffirst_precision_layer

Double_t splitcal::ffirst_precision_layer
private

Definition at line 118 of file splitcal.h.

◆ fLength

Double_t splitcal::fLength
private

time

Definition at line 107 of file splitcal.h.

◆ fMom

TLorentzVector splitcal::fMom
private

position at entrance

Definition at line 105 of file splitcal.h.

◆ fnECALSamplings

Int_t splitcal::fnECALSamplings
private

Definition at line 112 of file splitcal.h.

◆ fnHCALSamplings

Int_t splitcal::fnHCALSamplings
private

Definition at line 112 of file splitcal.h.

◆ fNModulesInX

Int_t splitcal::fNModulesInX
private

Definition at line 119 of file splitcal.h.

◆ fNModulesInY

Int_t splitcal::fNModulesInY
private

Definition at line 119 of file splitcal.h.

◆ fNStripsPerModule

Int_t splitcal::fNStripsPerModule
private

Definition at line 120 of file splitcal.h.

◆ fnum_precision_layers

Double_t splitcal::fnum_precision_layers
private

Definition at line 118 of file splitcal.h.

◆ fPos

TLorentzVector splitcal::fPos
private

volume id

Definition at line 104 of file splitcal.h.

◆ fsecond_precision_layer

Double_t splitcal::fsecond_precision_layer
private

Definition at line 118 of file splitcal.h.

◆ fsplitcalPointCollection

TClonesArray* splitcal::fsplitcalPointCollection
private

container for data points

Definition at line 125 of file splitcal.h.

◆ fStripHalfLength

Double_t splitcal::fStripHalfLength
private

Definition at line 121 of file splitcal.h.

◆ fStripHalfWidth

Double_t splitcal::fStripHalfWidth
private

Definition at line 121 of file splitcal.h.

◆ fthird_precision_layer

Double_t splitcal::fthird_precision_layer
private

Definition at line 118 of file splitcal.h.

◆ fTime

Double_t splitcal::fTime
private

momentum at entrance

Definition at line 106 of file splitcal.h.

◆ fTrackID

Int_t splitcal::fTrackID
private

Track information to be stored until the track leaves the active volume.

Definition at line 102 of file splitcal.h.

◆ fVolumeID

Int_t splitcal::fVolumeID
private

track index

Definition at line 103 of file splitcal.h.

◆ fXMax

Double_t splitcal::fXMax
private

Definition at line 116 of file splitcal.h.

◆ fYMax

Double_t splitcal::fYMax
private

Definition at line 117 of file splitcal.h.

◆ fZStart

Double_t splitcal::fZStart
private

Definition at line 114 of file splitcal.h.

◆ xfFilterECALThickness

Double_t splitcal::xfFilterECALThickness
private

Definition at line 109 of file splitcal.h.


The documentation for this class was generated from the following files: