PamVMC/src/PamVMCApplication.cxx

// $Id: PamVMCApplication.cxx,v 1.0 2007/06/01 
//
// Geant4 PAMELA VMC 
#include <cstdlib>
#include <string>
#include <iostream>
#include <TROOT.h>
#include <TInterpreter.h>
#include <TVirtualMC.h>
#include <Riostream.h>
#include <TGeoManager.h>
#include <TVirtualGeoTrack.h>
#include <TLorentzVector.h>
#include "PamVMCApplication.h"
#include "PamVMCStack.h"
#include "PamVMCPrimaryGenerator.h"

#include "PamVMCFieldMgr.h"

#include "PamVMCSDMgr.h"

#include "PamVMCVolCrossMgr.h"

#include "PamVMCDetPamela.h"

//#include "TFluka.h"


ClassImp(PamVMCApplication)

PamVMCApplication::PamVMCApplication(const char *name, const char *title, const char* filename, Int_t seed) 
  : TVirtualMCApplication(name,title),
    fEventNo(0),
    fVerbose(2),
    fStack(0),
    fdetector(0),
    fPrimaryGenerator(0),
    fRootManager(filename,kWrite)
{
// Standard constructor

  // Create a user stack
  fStack = new PamVMCStack(1000);
    
  // Create a primary generator
  fPrimaryGenerator = new PamVMCPrimaryGenerator(fStack);

  //Create random object
  frandom = new TRandom3(seed);

  fdetector = new PamVMCDetPamela();

  PamVMCFieldMgr::Instance()->LoadField();
  
  TString fname=filename; fname+=".pam";

  PamVMCRawMgr::Instance()-> CreateOutFile(fname.Data());
}


PamVMCApplication::~PamVMCApplication() 
{
// Destructor  
  delete frandom;
  delete fStack;
  delete fdetector;
  delete fPrimaryGenerator;
  delete gMC;
  gMC = 0;
}


void PamVMCApplication::InitMC(const char* setup)
{    
// Initialize MC.

  fVerbose.InitMC();

  gROOT->LoadMacro(setup);
  gInterpreter->ProcessLine("Config()");
  gMC->SetStack(fStack);
  gMC->SetRandom(frandom);
  
  std::cout<<"!!! BEFORE INIT GMC"<<std::endl;
  
  gMC->Init();
  cout<<"INIT MC"<<endl;
  gMC->BuildPhysics();
  cout<<"PHYSICS BUILDED"<<endl;

#if ROOT_VERSION_CODE >= 333572    
  gMC->SetMagField(PamVMCFieldMgr::Instance());
  cout<<"Magnetic field setted"<<endl;
#endif

  fdetector->InitMC();
  
  // fdetector->Print();
  
  PamVMCSDMgr::Instance()->Register();
  
  PamVMCDigMgr::Instance()->Initialize(gMC->GetRandom());

  PamVMCRawMgr::Instance()->WriteToFile();
  
  fPrimaryGenerator->SetRandom(gMC->GetRandom());
  fPrimaryGenerator->Register();
  
  std::cout<<"LIST OF ACTIVE SD'S:"<<std::endl;
  
  PamVMCSDMgr::Instance()->Print();
}

void PamVMCApplication::RunMC(Int_t nofEvents)
{    
   // MC run.

  fVerbose.RunMC(nofEvents);

  gMC->ProcessRun(nofEvents);

  fVerbose.FinishRun();
}

void PamVMCApplication::AddIons(){

  fVerbose.AddIons();
  //Arguments are:
  //name
  //Z (atomic number)
  //A (atomic mass)
  //Q (charge eplus)
  //exitation exitation energy [GeV] 0. if in ground state
  //mass [GeV], if not specified, calculated approx

  gMC->DefineIon("He3",2,3,2,0.,2.808);
  gMC->DefineIon("He4",2,4,2,0.,3.727); 

  //gMC->DefineIon("He7",2,7,2,0.);
  //gMC->DefineIon("Boron",5,10,5,0.);  
  //gMC->DefineIon("Carbon",6,12,6,0.);
}


void PamVMCApplication::ConstructGeometry()
{    
  // Construct geometry using detector contruction class

  fVerbose.ConstructGeometry();
  new TGeoManager("PAM_GEOM","PAMELA GEOMETRY");
  fdetector->ConstructGeometry();
  //TGeoManager::Import("/home/pamela/PamCAL210508/gp2tgeo.root");
  gGeoManager->CloseGeometry();
  gMC->SetRootGeometry();

  cout<<"!!CONSTRUCT GEOMETRY!!"<<endl;
  fdetector->ConstructCuts();

}

void PamVMCApplication::InitGeometry()
{    
  // Initialize geometry
  
  cout<<"!!SETUP CUTS"<<endl;
  //fdetector->InitGeometry();

  fVerbose.InitGeometry(); 
}


void PamVMCApplication::GeneratePrimaries()
{    
  // Fill the user stack (derived from TVirtualMCStack) with primary particles.
  
  fPrimaryGenerator->GeneratePrimary();
}

void PamVMCApplication::BeginEvent()
{    
// User actions at beginning of event
  PamVMCSDMgr::Instance()->PreEvent();
  fEventNo++;
}

void PamVMCApplication::BeginPrimary()
{    
// User actions at beginning of a primary track
  fPrimaryGenerator->SetGood(kTRUE);
  fVerbose.BeginPrimary();
}

void PamVMCApplication::PreTrack()
{    
// User actions at beginning of each track
 
}


#define CAVITY_X 8.07
#define CAVITY_Y 6.57
#define ZCAVITY_MIN 27.399
#define ZCAVITY_MAX 71.059


Bool_t PamVMCApplication::IsInsideCavity(){
  Double_t x,y,z;
  gMC->TrackPosition(x,y,z);
  if( (z<ZCAVITY_MIN) || (z>ZCAVITY_MAX)) return kTRUE;
  if( (TMath::Abs(x) > CAVITY_X) || (TMath::Abs(y) > CAVITY_Y)) return kFALSE;
  return kTRUE;
}


void PamVMCApplication::Stepping()
{    
  if (fStack->GetCurrentTrack()->IsPrimary()){
    if(!IsInsideCavity()) fPrimaryGenerator->SetGood(kFALSE);
  }
  
// User actions at each step
  fVerbose.Stepping();
  PamVMCDetectorSD* temp=(PamVMCDetectorSD*) PamVMCSDMgr::Instance()->GetSD(gMC->CurrentVolName());
  if(temp){
    fdstatus=INSIDE;
    if(gMC->IsTrackEntering()){ fdstatus=ENTERING; }
    if(gMC->IsTrackExiting() || gMC->IsTrackOut() || gMC->IsTrackStop() ) { fdstatus=EXITING;}
    temp->FillHit(fdstatus,gMC,fStack->GetCurrentTrack()->IsPrimary());    
  }
  
}

void PamVMCApplication::PostTrack()
{    
// User actions after finishing of each track

}

void PamVMCApplication::FinishPrimary()
{    
// User actions after finishing of a primary track
  //fPrimaryGenerator->SetGood(PamVMCVolCrossMgr::Instance()->IsTrackGood());
  //PamVMCVolCrossMgr::Instance()->Reset();
  fVerbose.FinishPrimary();
}

void PamVMCApplication::FinishEvent()
{    
// User actions after finishing of an event

    
    PamVMCSDMgr::Instance()->Digitize();
    //EventNo needs for trigger, particle PDG needs for TOF..
    PamVMCDigMgr::Instance()->Digitize(fEventNo,fPrimaryGenerator->GetParticle());

    PamVMCSDMgr::Instance()->Compress();

    PamVMCRawMgr::Instance()->WriteEvent();

    fRootManager.Fill();

    fPrimaryGenerator->ClearPrimCol();

    PamVMCSDMgr::Instance()->Clear();

    fStack->Reset();
   
} 

void PamVMCApplication::FinishRun()
{
  fVerbose.FinishRun();
  
  fRootManager.WriteAll();

  PamVMCDigMgr::Instance()->FinishRun(); //write RunTrailer to buffer

  PamVMCRawMgr::Instance()->FinishRun(); //write RunTrailer to file and close;

  PamVMCDigMgr::Instance()->PrintCollections();

}

#if ROOT_VERSION_CODE < 333572  
void PamVMCApplication::Field(const Double_t* x, Double_t* b) const
{
  PamVMCFieldMgr::Instance()->Field(x,b);
}
#endif
1	// $Id: PamVMCApplication.cxx,v 1.0 2007/06/01
2	//
3	// Geant4 PAMELA VMC
4	#include <cstdlib>
5	#include <string>
6	#include <iostream>
7	#include <TROOT.h>
8	#include <TInterpreter.h>
9	#include <TVirtualMC.h>
10	#include <Riostream.h>
11	#include <TGeoManager.h>
12	#include <TVirtualGeoTrack.h>
13	#include <TLorentzVector.h>
14	#include "PamVMCApplication.h"
15	#include "PamVMCStack.h"
16	#include "PamVMCPrimaryGenerator.h"
17
18	#include "PamVMCFieldMgr.h"
19
20	#include "PamVMCSDMgr.h"
21
22	#include "PamVMCVolCrossMgr.h"
23
24	#include "PamVMCDetPamela.h"
25
26	//#include "TFluka.h"
27
28
29
30	ClassImp(PamVMCApplication)
31
32	PamVMCApplication::PamVMCApplication(const char name, const char title, const char* filename, Int_t seed)
33	: TVirtualMCApplication(name,title),
34	fEventNo(0),
35	fVerbose(2),
36	fStack(0),
37	fdetector(0),
38	fPrimaryGenerator(0),
39	fRootManager(filename,kWrite)
40	{
41	// Standard constructor
42
43	// Create a user stack
44	fStack = new PamVMCStack(1000);
45
46	// Create a primary generator
47	fPrimaryGenerator = new PamVMCPrimaryGenerator(fStack);
48
49	//Create random object
50	frandom = new TRandom3(seed);
51
52	fdetector = new PamVMCDetPamela();
53
54	PamVMCFieldMgr::Instance()->LoadField();
55
56	TString fname=filename; fname+=".pam";
57
58	PamVMCRawMgr::Instance()-> CreateOutFile(fname.Data());
59	}
60
61
62	PamVMCApplication::~PamVMCApplication()
63	{
64	// Destructor
65	delete frandom;
66	delete fStack;
67	delete fdetector;
68	delete fPrimaryGenerator;
69	delete gMC;
70	gMC = 0;
71	}
72
73
74	void PamVMCApplication::InitMC(const char* setup)
75	{
76	// Initialize MC.
77
78	fVerbose.InitMC();
79
80	gROOT->LoadMacro(setup);
81	gInterpreter->ProcessLine("Config()");
82	gMC->SetStack(fStack);
83	gMC->SetRandom(frandom);
84
85	std::cout<<"!!! BEFORE INIT GMC"<<std::endl;
86
87	gMC->Init();
88	cout<<"INIT MC"<<endl;
89	gMC->BuildPhysics();
90	cout<<"PHYSICS BUILDED"<<endl;
91
92	#if ROOT_VERSION_CODE >= 333572
93	gMC->SetMagField(PamVMCFieldMgr::Instance());
94	cout<<"Magnetic field setted"<<endl;
95	#endif
96
97	fdetector->InitMC();
98
99	// fdetector->Print();
100
101	PamVMCSDMgr::Instance()->Register();
102
103	PamVMCDigMgr::Instance()->Initialize(gMC->GetRandom());
104
105	PamVMCRawMgr::Instance()->WriteToFile();
106
107	fPrimaryGenerator->SetRandom(gMC->GetRandom());
108	fPrimaryGenerator->Register();
109
110	std::cout<<"LIST OF ACTIVE SD'S:"<<std::endl;
111
112	PamVMCSDMgr::Instance()->Print();
113	}
114
115	void PamVMCApplication::RunMC(Int_t nofEvents)
116	{
117	// MC run.
118
119	fVerbose.RunMC(nofEvents);
120
121	gMC->ProcessRun(nofEvents);
122
123	fVerbose.FinishRun();
124	}
125
126	void PamVMCApplication::AddIons(){
127
128	fVerbose.AddIons();
129	//Arguments are:
130	//name
131	//Z (atomic number)
132	//A (atomic mass)
133	//Q (charge eplus)
134	//exitation exitation energy [GeV] 0. if in ground state
135	//mass [GeV], if not specified, calculated approx
136
137	gMC->DefineIon("He3",2,3,2,0.,2.808);
138	gMC->DefineIon("He4",2,4,2,0.,3.727);
139
140	//gMC->DefineIon("He7",2,7,2,0.);
141	//gMC->DefineIon("Boron",5,10,5,0.);
142	//gMC->DefineIon("Carbon",6,12,6,0.);
143	}
144
145
146	void PamVMCApplication::ConstructGeometry()
147	{
148	// Construct geometry using detector contruction class
149
150	fVerbose.ConstructGeometry();
151	new TGeoManager("PAM_GEOM","PAMELA GEOMETRY");
152	fdetector->ConstructGeometry();
153	//TGeoManager::Import("/home/pamela/PamCAL210508/gp2tgeo.root");
154	gGeoManager->CloseGeometry();
155	gMC->SetRootGeometry();
156
157	cout<<"!!CONSTRUCT GEOMETRY!!"<<endl;
158	fdetector->ConstructCuts();
159
160	}
161
162	void PamVMCApplication::InitGeometry()
163	{
164	// Initialize geometry
165
166	cout<<"!!SETUP CUTS"<<endl;
167	//fdetector->InitGeometry();
168
169	fVerbose.InitGeometry();
170	}
171
172
173	void PamVMCApplication::GeneratePrimaries()
174	{
175	// Fill the user stack (derived from TVirtualMCStack) with primary particles.
176
177	fPrimaryGenerator->GeneratePrimary();
178	}
179
180	void PamVMCApplication::BeginEvent()
181	{
182	// User actions at beginning of event
183	PamVMCSDMgr::Instance()->PreEvent();
184	fEventNo++;
185	}
186
187	void PamVMCApplication::BeginPrimary()
188	{
189	// User actions at beginning of a primary track
190	fPrimaryGenerator->SetGood(kTRUE);
191	fVerbose.BeginPrimary();
192	}
193
194	void PamVMCApplication::PreTrack()
195	{
196	// User actions at beginning of each track
197
198	}
199
200
201	#define CAVITY_X 8.07
202	#define CAVITY_Y 6.57
203	#define ZCAVITY_MIN 27.399
204	#define ZCAVITY_MAX 71.059
205
206
207	Bool_t PamVMCApplication::IsInsideCavity(){
208	Double_t x,y,z;
209	gMC->TrackPosition(x,y,z);
210	if( (z<ZCAVITY_MIN) \|\| (z>ZCAVITY_MAX)) return kTRUE;
211	if( (TMath::Abs(x) > CAVITY_X) \|\| (TMath::Abs(y) > CAVITY_Y)) return kFALSE;
212	return kTRUE;
213	}
214
215
216	void PamVMCApplication::Stepping()
217	{
218	if (fStack->GetCurrentTrack()->IsPrimary()){
219	if(!IsInsideCavity()) fPrimaryGenerator->SetGood(kFALSE);
220	}
221
222	// User actions at each step
223	fVerbose.Stepping();
224	PamVMCDetectorSD* temp=(PamVMCDetectorSD*) PamVMCSDMgr::Instance()->GetSD(gMC->CurrentVolName());
225	if(temp){
226	fdstatus=INSIDE;
227	if(gMC->IsTrackEntering()){ fdstatus=ENTERING; }
228	if(gMC->IsTrackExiting() \|\| gMC->IsTrackOut() \|\| gMC->IsTrackStop() ) { fdstatus=EXITING;}
229	temp->FillHit(fdstatus,gMC,fStack->GetCurrentTrack()->IsPrimary());
230	}
231
232	}
233
234	void PamVMCApplication::PostTrack()
235	{
236	// User actions after finishing of each track
237
238	}
239
240	void PamVMCApplication::FinishPrimary()
241	{
242	// User actions after finishing of a primary track
243	//fPrimaryGenerator->SetGood(PamVMCVolCrossMgr::Instance()->IsTrackGood());
244	//PamVMCVolCrossMgr::Instance()->Reset();
245	fVerbose.FinishPrimary();
246	}
247
248	void PamVMCApplication::FinishEvent()
249	{
250	// User actions after finishing of an event
251
252
253
254	PamVMCSDMgr::Instance()->Digitize();
255	//EventNo needs for trigger, particle PDG needs for TOF..
256	PamVMCDigMgr::Instance()->Digitize(fEventNo,fPrimaryGenerator->GetParticle());
257
258	PamVMCSDMgr::Instance()->Compress();
259
260	PamVMCRawMgr::Instance()->WriteEvent();
261
262	fRootManager.Fill();
263
264	fPrimaryGenerator->ClearPrimCol();
265
266	PamVMCSDMgr::Instance()->Clear();
267
268	fStack->Reset();
269
270	}
271
272	void PamVMCApplication::FinishRun()
273	{
274	fVerbose.FinishRun();
275
276	fRootManager.WriteAll();
277
278	PamVMCDigMgr::Instance()->FinishRun(); //write RunTrailer to buffer
279
280	PamVMCRawMgr::Instance()->FinishRun(); //write RunTrailer to file and close;
281
282	PamVMCDigMgr::Instance()->PrintCollections();
283
284	}
285
286	#if ROOT_VERSION_CODE < 333572
287	void PamVMCApplication::Field(const Double_t* x, Double_t* b) const
288	{
289	PamVMCFieldMgr::Instance()->Field(x,b);
290	}
291	#endif