PamVMC/examples/run_fluka_GP.C

// $Id: run_g4.C,v 1.0 2007/06/01 
//
// Macro for running G4PAM_VMC with 
class TFluka; 

void run_fluka_GP(const TString& configMacro = "flukaConfig.C") 
{
  // Load basic libraries
  gROOT->LoadMacro("../macro/basiclibs.C");
  basiclibs();

  // Load PAMELA libs
  TString PAMLIB=gSystem->Getenv("PAM_LIB");
  gSystem->Load(PAMLIB+"/libyoda.so");
  gSystem->Load(PAMLIB+"/libDarthVader.so");
  gSystem->Load(PAMLIB+"/libPamLevel2.so");


  // Load TFluka libraries
  gSystem->Load("/home/nikolas/francesco/vmc/fluka_vmc/lib/tgt_linux/libfluka");  
  // Load this example libraries
  TString G4WD=gSystem->Getenv("G4WORKDIR");
  TString PLAT=gSystem->Getenv("PLATFORM");
  gSystem->Load(G4WD+"/lib/tgt_"+PLAT+"/libPamVMC_fc");

  // MC application
  PamVMCApplication* appl 
    = new PamVMCApplication("PAMFLUKA_VMC", "PAMELA GEANT4 VMC application");
    

    appl->SetVerboseLevel(1);  
 
  
  // Initialize MC
  
  appl->InitMC(configMacro);
  // appl->SetVerboseLevel(1);  
  
  // Run MC
  

  //START READING INPUT FILE

  Int_t Ipa;
  Float_t X0, Y0, Z0, Theta, Phi, P0;
  TTree * h20 = 0;

  TString Filename = "pbar_r3_1-2.Beam.root";

  TString HOME=gSystem->Getenv("HOME");

  TSystemDirectory *readdir=new TSystemDirectory("readdir",HOME+"/RES");
    
  if ((TList*)(readdir->GetListOfFiles())->FindObject(Filename) ){
    cout<< Filename<<" found..."<<endl;
    TFile* inpf = new TFile(HOME+"/RES/"+Filename);
    h20 = (TTree*)inpf->Get("h20");
    if (h20->GetBranch("Ipa")){
      h20->SetBranchAddress("Ipa",&Ipa);
      h20->SetBranchStatus("Ipa",1);
    }
    if (h20->GetBranch("X0")){
      h20->SetBranchAddress("X0",&X0);
      h20->SetBranchStatus("X0",1);
    }
    if (h20->GetBranch("Y0")){
      h20->SetBranchAddress("Y0",&Y0);
      h20->SetBranchStatus("Y0",1);
    }
    if (h20->GetBranch("Z0")){
      h20->SetBranchAddress("Z0",&Z0);
      h20->SetBranchStatus("Z0",1);
    }
    if (h20->GetBranch("P0")){
      h20->SetBranchAddress("P0",&P0);
      h20->SetBranchStatus("P0",1);
    }
    if (h20->GetBranch("Theta")){
      h20->SetBranchAddress("Theta",&Theta);
      h20->SetBranchStatus("Theta",1);
    }
    if (h20->GetBranch("Phi")){
      h20->SetBranchAddress("Phi",&Phi);
      h20->SetBranchStatus("Phi",1);
    }
  }
  Int_t nevents = (Int_t)h20->GetEntries();

  //  nevents = 15;
  cout<<"PROCESS "<<nevents<<" EVENTS"<<endl;

  Double_t PX, PY, PZ;
  Int_t runNo = 0;
  for (Int_t i=0; i<nevents; i++){
    h20->GetEntry(i);
    cout<<"Ipa="<<Ipa<<endl;
    appl->GetPrimaryGenerator()->SetParticle(kProtonBar);
    appl->GetPrimaryGenerator()->SetMomentum(P0);
    cout<<"EKIN "<<appl->GetPrimaryGenerator()->GetKinEnergy()<<" GeV"<<endl;
    //cout<<"applP0="<<appl->GetPrimaryGenerator()->GetMomentum()<<" GV"<<endl;
    cout<<"applRIG "<<appl->GetPrimaryGenerator()->GetRigidity()<<" GV"<<endl
    PX=P0*sin(Theta)*cos(Phi);
    PY=P0*sin(Theta)*sin(Phi);
    PZ=-P0*cos(Theta);
    cout<<"(PX0,PY0,PZ0)="<<PX<<","<<PY<<","<<PZ<<endl;
    appl->GetPrimaryGenerator()->SetDirection(Theta,Phi);
    cout<<"(X0,Y0,Z0)"<<X0<<","<<Y0<<","<<Z0<<endl;
    appl->GetPrimaryGenerator()->SetPosition(X0,Y0,Z0);
    
    cout<<">StartRUN: "<<runNo<<endl;
    runNo++;

    appl->RunMC(1);

  }


  appl->FinishRun();
    delete appl;
  
  // return 
  
}  
1	// $Id: run_g4.C,v 1.0 2007/06/01
2	//
3	// Macro for running G4PAM_VMC with
4	class TFluka;
5
6	void run_fluka_GP(const TString& configMacro = "flukaConfig.C")
7	{
8	// Load basic libraries
9	gROOT->LoadMacro("../macro/basiclibs.C");
10	basiclibs();
11
12	// Load PAMELA libs
13	TString PAMLIB=gSystem->Getenv("PAM_LIB");
14	gSystem->Load(PAMLIB+"/libyoda.so");
15	gSystem->Load(PAMLIB+"/libDarthVader.so");
16	gSystem->Load(PAMLIB+"/libPamLevel2.so");
17
18
19	// Load TFluka libraries
20	gSystem->Load("/home/nikolas/francesco/vmc/fluka_vmc/lib/tgt_linux/libfluka");
21	// Load this example libraries
22	TString G4WD=gSystem->Getenv("G4WORKDIR");
23	TString PLAT=gSystem->Getenv("PLATFORM");
24	gSystem->Load(G4WD+"/lib/tgt_"+PLAT+"/libPamVMC_fc");
25
26	// MC application
27	PamVMCApplication* appl
28	= new PamVMCApplication("PAMFLUKA_VMC", "PAMELA GEANT4 VMC application");
29
30
31	appl->SetVerboseLevel(1);
32
33
34	// Initialize MC
35
36	appl->InitMC(configMacro);
37	// appl->SetVerboseLevel(1);
38
39	// Run MC
40
41
42	//START READING INPUT FILE
43
44	Int_t Ipa;
45	Float_t X0, Y0, Z0, Theta, Phi, P0;
46	TTree * h20 = 0;
47
48	TString Filename = "pbar_r3_1-2.Beam.root";
49
50	TString HOME=gSystem->Getenv("HOME");
51
52	TSystemDirectory *readdir=new TSystemDirectory("readdir",HOME+"/RES");
53
54	if ((TList*)(readdir->GetListOfFiles())->FindObject(Filename) ){
55	cout<< Filename<<" found..."<<endl;
56	TFile* inpf = new TFile(HOME+"/RES/"+Filename);
57	h20 = (TTree*)inpf->Get("h20");
58	if (h20->GetBranch("Ipa")){
59	h20->SetBranchAddress("Ipa",&Ipa);
60	h20->SetBranchStatus("Ipa",1);
61	}
62	if (h20->GetBranch("X0")){
63	h20->SetBranchAddress("X0",&X0);
64	h20->SetBranchStatus("X0",1);
65	}
66	if (h20->GetBranch("Y0")){
67	h20->SetBranchAddress("Y0",&Y0);
68	h20->SetBranchStatus("Y0",1);
69	}
70	if (h20->GetBranch("Z0")){
71	h20->SetBranchAddress("Z0",&Z0);
72	h20->SetBranchStatus("Z0",1);
73	}
74	if (h20->GetBranch("P0")){
75	h20->SetBranchAddress("P0",&P0);
76	h20->SetBranchStatus("P0",1);
77	}
78	if (h20->GetBranch("Theta")){
79	h20->SetBranchAddress("Theta",&Theta);
80	h20->SetBranchStatus("Theta",1);
81	}
82	if (h20->GetBranch("Phi")){
83	h20->SetBranchAddress("Phi",&Phi);
84	h20->SetBranchStatus("Phi",1);
85	}
86	}
87	Int_t nevents = (Int_t)h20->GetEntries();
88
89	// nevents = 15;
90	cout<<"PROCESS "<<nevents<<" EVENTS"<<endl;
91
92	Double_t PX, PY, PZ;
93	Int_t runNo = 0;
94	for (Int_t i=0; i<nevents; i++){
95	h20->GetEntry(i);
96	cout<<"Ipa="<<Ipa<<endl;
97	appl->GetPrimaryGenerator()->SetParticle(kProtonBar);
98	appl->GetPrimaryGenerator()->SetMomentum(P0);
99	cout<<"EKIN "<<appl->GetPrimaryGenerator()->GetKinEnergy()<<" GeV"<<endl;
100	//cout<<"applP0="<<appl->GetPrimaryGenerator()->GetMomentum()<<" GV"<<endl;
101	cout<<"applRIG "<<appl->GetPrimaryGenerator()->GetRigidity()<<" GV"<<endl
102	PX=P0sin(Theta)cos(Phi);
103	PY=P0sin(Theta)sin(Phi);
104	PZ=-P0*cos(Theta);
105	cout<<"(PX0,PY0,PZ0)="<<PX<<","<<PY<<","<<PZ<<endl;
106	appl->GetPrimaryGenerator()->SetDirection(Theta,Phi);
107	cout<<"(X0,Y0,Z0)"<<X0<<","<<Y0<<","<<Z0<<endl;
108	appl->GetPrimaryGenerator()->SetPosition(X0,Y0,Z0);
109
110	cout<<">StartRUN: "<<runNo<<endl;
111	runNo++;
112
113	appl->RunMC(1);
114
115	}
116
117
118
119
120	appl->FinishRun();
121	delete appl;
122
123	// return
124
125	}