[u/mrichter/AliRoot.git] / ACORDE / macros / TestACORDEHits.C

/////////////////////////////////////////////////////////////////////////
//   This macro is a small example of a ROOT macro
//   illustrating how to read the hits in ACORDE,
//   fill some histograms and store it on a file
//
//   Author:
//
//	Eleazar Cuautle <ecuautle@nucleares.unam.mx>
//
//   Also:
//
//	Mario Rodriguez Cahuantzi <mrodrigu@mail.cern.ch>
//             
/////////////////////////////////////////////////////////////////////////

void TestACORDEHits (const char *filename="galice.root",Int_t numberOfEvents=0){
   // delete existing gAlice
  // Dynamically link some shared libs
  if (gClassTable->GetID("AliRun") < 0) 
    {
      gROOT->LoadMacro("loadlibs.C");
      loadlibs();
    }
  if (gAlice)
    {
      //delete gAlice->GetRunLoader();
      delete gAlice;
      gAlice = 0x0;
    }
  
  AliRunLoader *rl = AliRunLoader::Open("galice.root",
					AliConfig::GetDefaultEventFolderName(),"read");
  if (!rl)
    {
      cerr<<"Can't load RunLoader from file! \n";      return 0x0;
    }
  
  rl->LoadgAlice();
  gAlice = rl->GetAliRun();
  if (!gAlice)
    {
      cerr << " AliRun object not found on file \n";     return 0x0; 
    }
  rl->LoadHeader();
  
  // Get the pointer to the ACORDE detector
  AliLoader *acordel = rl->GetLoader("ACORDELoader");
  AliACORDE *ACORDE = (AliACORDE *) gAlice->GetDetector("ACORDE");
  if (ACORDE == 0x0 || acordel == 0x0) {
    cerr << " Can not find ACORDE or ACORDELoader \n";    return 0x0;  
  }
  //Determine how many events would you like to analyze
  if (numberOfEvents==0) numberOfEvents=(Int_t)(rl->GetNumberOfEvents());
  
  //create histograms

  TH1F	*fACORDEEnergy = new TH1F("fACORDEEnergy","Energy distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Energy (GeV)",1200,0,360);
  TH1F	*fACORDEEnergyLoss = new TH1F("fACORDEEnergyLoss","Loss energy distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);Energy loss (Gev)",100,0,0.030);
  TH1F	*fACORDEAzimuth = new TH1F("fACORDEAzimuth","Azimuth angle distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Azimuth angle",360,-180,180);
  TH1F	*fACORDEPolar = new TH1F("fACORDEPolar","Polar angle distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Polar angle",180,0,180);
  TH2F  *fACORDExz = new TH2F("fACORDExz" ,"Distribution XZ  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); X (cms); Z (cms) ",900,-850.,850.,1200,-500.,500.); 
  TH2F  *fACORDEAzimPol =  new TH2F("fACORDEAzimPol" ,"Azimuth vs Polar angle distribution ",360,-180.,180.,180,0.,180.);
  TH1F  *fACORDEPx = new TH1F("fACORDEPx" ,"P_{x} distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{x}",60,-120.,120.);
  TH1F  *fACORDEPy = new TH1F("fACORDEPy" ,"P_{y} distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{y}",60,-120.,120.);
  TH1F  *fACORDEPz = new TH1F("fACORDEPz" ,"P_{z} distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{z}",60,-120.,120.);
  TH1F  *fACORDEPt = new TH1F("fACORDEPt" ,"P_{t} distribution  of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{t}",1200,0.,360.);
  TH2F  *fELossAzimuth = new TH2F("fELossAzimuth","Dist. Azimuth Angle VS Energy Loss; Azimuth angle; Energy (GeV)",360,0,350,100,0,0.030); 
  TH2F  *fEnergyAzimuth = new TH2F("fEnergyAzimuth","Dist. Azimuth Angle VS Energy; Azimuth angle; Energy (GeV)",360,0,360,300,0,350); 
  TH2F  *fELossPolar = new TH2F("fELossPolar","Dist. Polar Angle VS Energy Loss; Polar angle; Energy (GeV)",180,0,180,100,0,0.030); 
  TH2F  *fEnergyPolar = new TH2F("fEnergyPolar","Dist. Polar Angle VS Energy; Polar angle; Energy (GeV)",180,0,180,300,0,350); 

 //------------------------------------------------------------------//
  
  for (Int_t ievent=0; ievent<numberOfEvents; ievent++) {
    if ((ievent%10) == 0)  printf ("Processing event %d \n", ievent);
    rl->GetEvent(ievent);
    
    // Get the pointer Hit tree
    acordel->LoadHits();
    TTree *hitTree = acordel->TreeH();
    ACORDE->SetTreeAddress();
    if (!hitTree) {
      cout << " No TreeH found" << endl;      return 0x0; //rc;
    }
    
    rl->LoadKinematics();
    Int_t nTrack = (Int_t) hitTree->GetEntries();
    //    cout <<"<AliACORDEanalyzeHits> Found "<< nTrack <<" primary particles with hits \n";
    
    // Start loop on tracks in the hits containers
    for(Int_t iTrack=0; iTrack<nTrack;iTrack++){
      ACORDE->ResetHits();
      hitTree->GetEvent(iTrack);          
      if(ACORDE) {	 
	for(acordeHit=(AliACORDEhit*)ACORDE->FirstHit(-1);acordeHit;
	    acordeHit=(AliACORDEhit*)ACORDE->NextHit()) {
		fACORDEEnergy->Fill(acordeHit->Energy());
		fACORDEEnergyLoss->Fill(acordeHit->Eloss());
		fACORDEAzimuth->Fill(acordeHit->AzimuthAngle());
		fACORDEPolar->Fill(acordeHit->PolarAngle());
		fACORDExz->Fill(acordeHit->X(),acordeHit->Z());
         	fELossAzimuth->Fill((Float_t)(acordeHit->AzimuthAngle()),(Float_t)(acordeHit->Eloss()));
 		fACORDEPx->Fill( (Float_t)(acordeHit->Px()));
	 	fACORDEPy->Fill( (Float_t)(acordeHit->Py()));
	  	fACORDEPz->Fill( (Float_t)(acordeHit->Pz()));
	 	fACORDEPt->Fill( TMath::Sqrt( (acordeHit->Px())*(acordeHit->Px())+
					     (acordeHit->Py())*(acordeHit->Py())) );
	  	fACORDEAzimPol->Fill( (Float_t)(acordeHit->AzimuthAngle()),
			        (Float_t)(acordeHit->PolarAngle()));
		fEnergyAzimuth->Fill(acordeHit->Energy(),acordeHit->AzimuthAngle());
		fEnergyPolar->Fill(acordeHit->Energy(),acordeHit->PolarAngle());
		fELossPolar->Fill(acordeHit->PolarAngle(),acordeHit->Eloss());
	//cout << "Muon triggered by ACORDE "<<endl;
	//cout << "Energy: "<<acordeHit->Energy()<<" EnergyLoss: "<<acordeHit->Eloss()<<" Polar angle: "<<acordeHit->PolarAngle()<<" Azimuth angle: "<<acordeHit->AzimuthAngle()<<endl;

	} // end for acordeHits
      } // end if ACORDE
    } // end for iTrack  
  }

 // Put the histograms in a TCanvas

	TCanvas *c1 = new TCanvas("c1","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
	c1->Divide(2,1);

	c1->cd(1)->SetLogy();
	fACORDEEnergy->SetMarkerStyle(kFullCircle);
	fACORDEEnergy->SetMarkerSize(0.7);
	fACORDEEnergy->Draw("E");

	c1->cd(2)->SetLogy();
	fACORDEEnergyLoss->SetMarkerStyle(kFullCircle);
	fACORDEEnergyLoss->SetMarkerSize(0.7);
	fACORDEEnergyLoss->Draw("E");
	
	c1->Draw();

	TCanvas *c2 = new TCanvas("c2","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
	c2->Divide(2,2);
	c2->cd(1);
	fACORDEAzimuth->Draw();

	c2->cd(2);
	fACORDEPolar->Draw();

	c2->cd(3);
	gStyle->SetPalette(1);
	fACORDExz->DrawCopy("colz");
	
	c2->cd(4);
	gStyle->SetPalette(1);
	fACORDEAzimPol->DrawCopy("colz");

	c2->Draw();

	TCanvas *c3 = new TCanvas("c3","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
	c3->Divide(2,2);

	c3->cd(1)->SetLogy();
	fACORDEPx->SetMarkerStyle(kFullCircle);
	fACORDEPx->SetMarkerSize(0.7);
	fACORDEPx->Draw("E");

	c3->cd(2)->SetLogy();
	fACORDEPy->SetMarkerStyle(kFullCircle);
	fACORDEPy->SetMarkerSize(0.7);
	fACORDEPy->Draw("E");

	c3->cd(3)->SetLogy();
	fACORDEPz->SetMarkerStyle(kFullCircle);
	fACORDEPz->SetMarkerSize(0.7);
	fACORDEPz->Draw("E");

	c3->cd(4)->SetLogy();
	fACORDEPt->SetMarkerStyle(kFullCircle);
	fACORDEPt->SetMarkerSize(0.7);
	fACORDEPt->Draw("E");

	c3->Draw();

	TCanvas *c4 = new TCanvas("c4","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
	c4->Divide(2,2);

	c4->cd(1);
	gStyle->SetPalette(1);
	fEnergyPolar->DrawCopy("colz");

	c4->cd(2);
	gStyle->SetPalette(1);
	fEnergyAzimuth->DrawCopy("colz");

	c4->cd(3);
	gStyle->SetPalette(1);
	fELossPolar->DrawCopy("colz");

	c4->cd(4);
	gStyle->SetPalette(1);
	fELossAzimuth->DrawCopy("colz");
	
	c4->Draw();


 // save histos in a root file
  TFile *fout = new TFile("ACORDE_hits.root","RECREATE");
	fACORDEEnergy->Write();
	fACORDEEnergyLoss->Write();
	fACORDEAzimuth->Write();
	fACORDEPolar->Write();
	fACORDExz->Write();
	fELossAzimuth->Write();
	fACORDEPx->Write();
	fACORDEPy->Write();
	fACORDEPz->Write();
	fACORDEPt->Write();
	fACORDEAzimPol->Write();
	c1->Write();
	c2->Write();
	c3->Write();
	c4->Write();
}
Commit	Line	Data
b5f2577c	1	/////////////////////////////////////////////////////////////////////////
	2	// This macro is a small example of a ROOT macro
	3	// illustrating how to read the hits in ACORDE,
	4	// fill some histograms and store it on a file
	5	//
	6	// Author:
	7	//
	8	// Eleazar Cuautle <ecuautle@nucleares.unam.mx>
	9	//
	10	// Also:
	11	//
	12	// Mario Rodriguez Cahuantzi <mrodrigu@mail.cern.ch>
	13	//
	14	/////////////////////////////////////////////////////////////////////////
	15
	16	void TestACORDEHits (const char *filename="galice.root",Int_t numberOfEvents=0){
	17	// delete existing gAlice
	18	// Dynamically link some shared libs
	19	if (gClassTable->GetID("AliRun") < 0)
	20	{
	21	gROOT->LoadMacro("loadlibs.C");
	22	loadlibs();
	23	}
	24	if (gAlice)
	25	{
	26	//delete gAlice->GetRunLoader();
	27	delete gAlice;
	28	gAlice = 0x0;
	29	}
	30
	31	AliRunLoader *rl = AliRunLoader::Open("galice.root",
	32	AliConfig::GetDefaultEventFolderName(),"read");
	33	if (!rl)
	34	{
	35	cerr<<"Can't load RunLoader from file! \n"; return 0x0;
	36	}
	37
	38	rl->LoadgAlice();
	39	gAlice = rl->GetAliRun();
	40	if (!gAlice)
	41	{
	42	cerr << " AliRun object not found on file \n"; return 0x0;
	43	}
	44	rl->LoadHeader();
	45
	46	// Get the pointer to the ACORDE detector
	47	AliLoader *acordel = rl->GetLoader("ACORDELoader");
	48	AliACORDE ACORDE = (AliACORDE ) gAlice->GetDetector("ACORDE");
	49	if (ACORDE == 0x0 \|\| acordel == 0x0) {
	50	cerr << " Can not find ACORDE or ACORDELoader \n"; return 0x0;
	51	}
	52	//Determine how many events would you like to analyze
	53	if (numberOfEvents==0) numberOfEvents=(Int_t)(rl->GetNumberOfEvents());
	54
	55	//create histograms
	56
	57	TH1F *fACORDEEnergy = new TH1F("fACORDEEnergy","Energy distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Energy (GeV)",1200,0,360);
	58	TH1F *fACORDEEnergyLoss = new TH1F("fACORDEEnergyLoss","Loss energy distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);Energy loss (Gev)",100,0,0.030);
	59	TH1F *fACORDEAzimuth = new TH1F("fACORDEAzimuth","Azimuth angle distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Azimuth angle",360,-180,180);
	60	TH1F *fACORDEPolar = new TH1F("fACORDEPolar","Polar angle distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); Polar angle",180,0,180);
	61	TH2F *fACORDExz = new TH2F("fACORDExz" ,"Distribution XZ of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE); X (cms); Z (cms) ",900,-850.,850.,1200,-500.,500.);
	62	TH2F *fACORDEAzimPol = new TH2F("fACORDEAzimPol" ,"Azimuth vs Polar angle distribution ",360,-180.,180.,180,0.,180.);
	63	TH1F *fACORDEPx = new TH1F("fACORDEPx" ,"P_{x} distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{x}",60,-120.,120.);
	64	TH1F *fACORDEPy = new TH1F("fACORDEPy" ,"P_{y} distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{y}",60,-120.,120.);
65	TH1F *fACORDEPz = new TH1F("fACORDEPz" ,"P_{z} distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{z}",60,-120.,120.);
66	TH1F *fACORDEPt = new TH1F("fACORDEPt" ,"P_{t} distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE);GeV/c;dN/dP_{t}",1200,0.,360.);
67	TH2F *fELossAzimuth = new TH2F("fELossAzimuth","Dist. Azimuth Angle VS Energy Loss; Azimuth angle; Energy (GeV)",360,0,350,100,0,0.030);
68	TH2F *fEnergyAzimuth = new TH2F("fEnergyAzimuth","Dist. Azimuth Angle VS Energy; Azimuth angle; Energy (GeV)",360,0,360,300,0,350);
69	TH2F *fELossPolar = new TH2F("fELossPolar","Dist. Polar Angle VS Energy Loss; Polar angle; Energy (GeV)",180,0,180,100,0,0.030);
70	TH2F *fEnergyPolar = new TH2F("fEnergyPolar","Dist. Polar Angle VS Energy; Polar angle; Energy (GeV)",180,0,180,300,0,350);
71
72	//------------------------------------------------------------------//
73
74	for (Int_t ievent=0; ievent<numberOfEvents; ievent++) {
75	if ((ievent%10) == 0) printf ("Processing event %d \n", ievent);
76	rl->GetEvent(ievent);
77
78	// Get the pointer Hit tree
79	acordel->LoadHits();
80	TTree *hitTree = acordel->TreeH();
81	ACORDE->SetTreeAddress();
82	if (!hitTree) {
83	cout << " No TreeH found" << endl; return 0x0; //rc;
84	}
85
86	rl->LoadKinematics();
87	Int_t nTrack = (Int_t) hitTree->GetEntries();
88	// cout <<"<AliACORDEanalyzeHits> Found "<< nTrack <<" primary particles with hits \n";
89
90	// Start loop on tracks in the hits containers
91	for(Int_t iTrack=0; iTrack<nTrack;iTrack++){
92	ACORDE->ResetHits();
93	hitTree->GetEvent(iTrack);
94	if(ACORDE) {
95	for(acordeHit=(AliACORDEhit*)ACORDE->FirstHit(-1);acordeHit;
96	acordeHit=(AliACORDEhit*)ACORDE->NextHit()) {
97	fACORDEEnergy->Fill(acordeHit->Energy());
98	fACORDEEnergyLoss->Fill(acordeHit->Eloss());
99	fACORDEAzimuth->Fill(acordeHit->AzimuthAngle());
100	fACORDEPolar->Fill(acordeHit->PolarAngle());
101	fACORDExz->Fill(acordeHit->X(),acordeHit->Z());
102	fELossAzimuth->Fill((Float_t)(acordeHit->AzimuthAngle()),(Float_t)(acordeHit->Eloss()));
103	fACORDEPx->Fill( (Float_t)(acordeHit->Px()));
104	fACORDEPy->Fill( (Float_t)(acordeHit->Py()));
105	fACORDEPz->Fill( (Float_t)(acordeHit->Pz()));
106	fACORDEPt->Fill( TMath::Sqrt( (acordeHit->Px())*(acordeHit->Px())+
107	(acordeHit->Py())*(acordeHit->Py())) );
108	fACORDEAzimPol->Fill( (Float_t)(acordeHit->AzimuthAngle()),
109	(Float_t)(acordeHit->PolarAngle()));
110	fEnergyAzimuth->Fill(acordeHit->Energy(),acordeHit->AzimuthAngle());
111	fEnergyPolar->Fill(acordeHit->Energy(),acordeHit->PolarAngle());
112	fELossPolar->Fill(acordeHit->PolarAngle(),acordeHit->Eloss());
113	//cout << "Muon triggered by ACORDE "<<endl;
114	//cout << "Energy: "<<acordeHit->Energy()<<" EnergyLoss: "<<acordeHit->Eloss()<<" Polar angle: "<<acordeHit->PolarAngle()<<" Azimuth angle: "<<acordeHit->AzimuthAngle()<<endl;
115
116	} // end for acordeHits
117	} // end if ACORDE
118	} // end for iTrack
119	}
120
121	// Put the histograms in a TCanvas
122
123	TCanvas *c1 = new TCanvas("c1","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
124	c1->Divide(2,1);
125
126	c1->cd(1)->SetLogy();
127	fACORDEEnergy->SetMarkerStyle(kFullCircle);
128	fACORDEEnergy->SetMarkerSize(0.7);
129	fACORDEEnergy->Draw("E");
130
131	c1->cd(2)->SetLogy();
132	fACORDEEnergyLoss->SetMarkerStyle(kFullCircle);
133	fACORDEEnergyLoss->SetMarkerSize(0.7);
134	fACORDEEnergyLoss->Draw("E");
135
136	c1->Draw();
137
138	TCanvas *c2 = new TCanvas("c2","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
139	c2->Divide(2,2);
140	c2->cd(1);
141	fACORDEAzimuth->Draw();
142
143	c2->cd(2);
144	fACORDEPolar->Draw();
145
146	c2->cd(3);
147	gStyle->SetPalette(1);
148	fACORDExz->DrawCopy("colz");
149
150	c2->cd(4);
151	gStyle->SetPalette(1);
152	fACORDEAzimPol->DrawCopy("colz");
153
154	c2->Draw();
155
156	TCanvas *c3 = new TCanvas("c3","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
157	c3->Divide(2,2);
158
159	c3->cd(1)->SetLogy();
160	fACORDEPx->SetMarkerStyle(kFullCircle);
161	fACORDEPx->SetMarkerSize(0.7);
162	fACORDEPx->Draw("E");
163
164	c3->cd(2)->SetLogy();
165	fACORDEPy->SetMarkerStyle(kFullCircle);
166	fACORDEPy->SetMarkerSize(0.7);
167	fACORDEPy->Draw("E");
168
169	c3->cd(3)->SetLogy();
170	fACORDEPz->SetMarkerStyle(kFullCircle);
171	fACORDEPz->SetMarkerSize(0.7);
172	fACORDEPz->Draw("E");
173
174	c3->cd(4)->SetLogy();
175	fACORDEPt->SetMarkerStyle(kFullCircle);
176	fACORDEPt->SetMarkerSize(0.7);
177	fACORDEPt->Draw("E");
178
179	c3->Draw();
180
181	TCanvas *c4 = new TCanvas("c4","ACORDE Hits distribution of muons generated by AliGenCosmicParam at z = 900 cms (Just above ACORDE)");
182	c4->Divide(2,2);
183
184	c4->cd(1);
185	gStyle->SetPalette(1);
186	fEnergyPolar->DrawCopy("colz");
187
188	c4->cd(2);
189	gStyle->SetPalette(1);
190	fEnergyAzimuth->DrawCopy("colz");
191
192	c4->cd(3);
193	gStyle->SetPalette(1);
194	fELossPolar->DrawCopy("colz");
195
196	c4->cd(4);
197	gStyle->SetPalette(1);
198	fELossAzimuth->DrawCopy("colz");
199
200	c4->Draw();
201
202
203	// save histos in a root file
204	TFile *fout = new TFile("ACORDE_hits.root","RECREATE");
205	fACORDEEnergy->Write();
206	fACORDEEnergyLoss->Write();
207	fACORDEAzimuth->Write();
208	fACORDEPolar->Write();
209	fACORDExz->Write();
210	fELossAzimuth->Write();
211	fACORDEPx->Write();
212	fACORDEPy->Write();
213	fACORDEPz->Write();
214	fACORDEPt->Write();
215	fACORDEAzimPol->Write();
216	c1->Write();
217	c2->Write();
218	c3->Write();
219	c4->Write();
220	}