TuneA settings removed from kPyJets
[u/mrichter/AliRoot.git] / EVGEN / TestPrimaries.C
CommitLineData
754a67f2 1void TestPrimaries(Int_t evNumber1=0, Int_t evNumber2=0)
2{
3// Dynamically link some shared libs
4 if (gClassTable->GetID("AliRun") < 0) {
5 gROOT->LoadMacro("loadlibs.C");
6 loadlibs();
7 }
8// Connect the Root Galice file containing Geometry, Kine and Hits
1d93d37b 9 AliRunLoader* rl = AliRunLoader::Open("galice.root");
754a67f2 10//
1d93d37b 11 TDatabasePDG* DataBase = TDatabasePDG::Instance();
754a67f2 12
13// Create some histograms
14
15 TH1F *thetaH = new TH1F("thetaH","Theta distribution",180,0,180);
ccbdd4de 16 thetaH->Sumw2();
754a67f2 17 TH1F *phiH = new TH1F("phiH","Phi distribution" ,180,-180,180);
ccbdd4de 18 phiH->Sumw2();
19 TH1F *phiH1 = new TH1F("phiH1","Phi distribution" ,180,-180,180);
20 phiH1->Sumw2();
754a67f2 21 TH1F *etaH = new TH1F("etaH","Pseudorapidity",120,-12,12);
ccbdd4de 22 etaH->Sumw2();
754a67f2 23 TH1F *yH = new TH1F("yH","Rapidity distribution",120,-12,12);
ccbdd4de 24 yH->Sumw2();
754a67f2 25 TH1F *eH = new TH1F("eH","Energy distribution",100,0,100);
ccbdd4de 26 eH->Sumw2();
a78c8623 27 TH1F *eetaH = new TH1F("eetaH","Pseudorapidity",120,0,12);
ccbdd4de 28 eetaH->Sumw2();
754a67f2 29 TH1F *ptH = new TH1F("ptH","Pt distribution",150,0,15);
ccbdd4de 30 ptH->Sumw2();
754a67f2 31//
32// Loop over events
33//
1d93d37b 34 rl->LoadKinematics();
35 rl->LoadHeader();
754a67f2 36
37 for (Int_t nev=0; nev<= evNumber2; nev++) {
1d93d37b 38 rl->GetEvent(nev);
39 AliStack* stack = rl->Stack();
40 Int_t npart = stack->GetNprimary();
754a67f2 41 if (nev < evNumber1) continue;
754a67f2 42//
43// Loop over primary particles (jpsi. upsilon, ...)
44//
14b783fa 45
754a67f2 46 for (Int_t part=0; part<npart; part++) {
1d93d37b 47 TParticle *MPart = stack->Particle(part);
754a67f2 48 Int_t mpart = MPart->GetPdgCode();
49 Int_t child1 = MPart->GetFirstDaughter();
50 Int_t child2 = MPart->GetLastDaughter();
51 Int_t mother = MPart->GetFirstMother();
754a67f2 52 Float_t Pt = MPart->Pt();
53 Float_t E = MPart->Energy();
54 Float_t Pz = MPart->Pz();
55 Float_t Py = MPart->Py();
56 Float_t Px = MPart->Px();
57 Float_t pT = TMath::Sqrt(Px*Px+Py*Py);
58 Float_t theta = MPart->Theta();
59 Float_t phi = MPart->Phi()-TMath::Pi();
60 Float_t eta = -TMath::Log(TMath::Tan(theta/2.));
73e54f72 61 Float_t y = 0.5*TMath::Log((E+Pz+1.e-13)/(E-Pz+1.e-13));
a78c8623 62
34e87c61 63 if (child1 >= 0) continue;
1d93d37b 64 if (DataBase->GetParticle(mpart)->Charge() == 0) continue;
a78c8623 65 Float_t wgt = 1./(Float_t ((evNumber2-evNumber1)+1.));
66 thetaH->Fill(theta*180./TMath::Pi(),wgt);
ccbdd4de 67 if (pT<2)
68 phiH->Fill(phi*180./TMath::Pi(),wgt);
69 else
70 phiH1->Fill(phi*180./TMath::Pi(),wgt);
a78c8623 71 etaH->Fill(eta,5.*wgt);
754a67f2 72 eetaH->Fill(eta,E);
a78c8623 73 yH->Fill(y,5.*wgt);
74 eH->Fill(E,wgt);
75 ptH->Fill(pT,wgt);
754a67f2 76 } // primary loop
77 }
78
79//Create a canvas, set the view range, show histograms
80 TCanvas *c1 = new TCanvas("c1","Canvas 1",400,10,600,700);
81 c1->Divide(2,2);
82 c1->cd(1); ptH->Draw();
83 c1->cd(2); etaH->Draw();
84 c1->cd(3); yH->Draw();
1d93d37b 85 c1->cd(4); eH->Draw();
754a67f2 86
87 TCanvas *c2 = new TCanvas("c2","Canvas 1",400,10,600,700);
88 c2->Divide(2,2);
89 c2->cd(1); phiH->Draw();
ccbdd4de 90 c2->cd(2); phiH1->Draw();
91 c2->cd(3); thetaH->Draw();
92 c2->cd(4); eetaH->Draw();
a78c8623 93
1d93d37b 94 TFile* f = new TFile("kineH.root", "update");
95 ptH->Write();
96 etaH->Write();
97 yH->Write();
98 phiH->Write();
ccbdd4de 99 phiH1->Write();
1d93d37b 100 f->Close();
101
a78c8623 102
754a67f2 103}
104
105
106
107
108
109
110