2 // This macro is a part of "Alice PPR fast flow analysis package"
5 // 1) open existing database with Flow Pico DST ("flowPicoEvent.root")
6 // 2) open source data - galice.root type file with events
7 // 3) reconstruct psi, psiA, psiB (sub-events in eta space) and v2.
8 // 4) fill Flow Pico DST with reconstructed data
11 // inputName: file with events
12 // if (inputName == 0) macro tries to fetch system variable FILE_NAME
13 // and use its value as the name of the file. (useful for bath processing)
15 // maxN: maximum number of events to be processed, if 0 process all events
18 // file "flowPicoEvent.root" with "flowData" NTuple have to exist
19 // it can be created by "AliFlowCreate.C"
21 // Sylwester Radomski, GSI
22 // mail: S.Radomski@gsi
25 ////////////////////////////////////////////////////////////////////////////////
28 //#include "TNtuple.h"
32 void AliFlowReconstruction(const char* inputName = 0, Int_t maxN = 0) {
34 const char *dataName = "flowPicoEvent.root";
39 Double_t trueV2, truePsi;
40 Double_t v2, psi, psiA, psiB;
44 if (gAlice) delete gAlice;
47 TFile *dataFile = new TFile(dataName,"UPDATE");
48 data = (TNtuple *) dataFile->Get("flowData");
53 TFile *inputFile = new TFile(gSystem->Getenv("FILE_NAME"), "READ");
54 ::Info("AliFlowReconstruction", "Using: %s",gSystem->Getenv("FILE_NAME"));
56 TFile *inputFile = new TFile(inputName, "READ");
59 gAlice = (AliRun*) inputFile->Get("gAlice");
61 N = gAlice->GetEventsPerRun();
62 if (maxN != 0 && maxN < N) N = maxN;
64 ::Info("AliFlowReconstruction", "Number of events to be processed = %d", N);
66 for (Int_t i=0; i<N; i++) {
68 ::Info("Processing event i = %d", i);
72 AliHeader *header = gAlice->GetHeader();
73 AliGenGeVSimEventHeader *egHeader = (AliGenGeVSimEventHeader*)header->GenEventHeader();
75 truePsi = egHeader->GetEventPlane() * 180 / TMath::Pi();
76 while(truePsi > 90) truePsi -= 180;
78 trueV2 = egHeader->GetEllipticFlow(211);
80 mult = gAlice->GetNtrack();
81 kine = gAlice->TreeK();
83 psi = FlowPsi(kine, 0);
84 psiA = FlowPsi(kine, 1);
85 psiB = FlowPsi(kine, -1);
87 v2 = FlowV2(kine, psi);
88 data->Fill(1, i, mult, truePsi, trueV2, psi, psiA, psiB, v2);
108 ////////////////////////////////////////////////////////////////////////////////
110 Double_t FlowPsi(TTree *kine, Int_t etaPart) {
112 // this function calculates psi from list of particles
113 // etaPart - division to sub-events in eta
119 Int_t mult = kine->GetEntries();
123 Double_t Ssin = 0, Scos = 0;
125 kine->Draw("Particles->Phi():Particles->Eta()", "", "goff");
129 for (Int_t i=0; i<mult; i++) {
131 if ( (etaPart * eta[i]) < 0 ) continue;
133 Ssin += sin( 2 * phi[i] );
134 Scos += cos( 2 * phi[i] );
137 //psi = atan( Ssin / Scos ) / 2;
138 psi = atan2 (Ssin, Scos) / 2;
139 psi = psi * 180. / TMath::Pi();
144 ////////////////////////////////////////////////////////////////////////////////
146 Double_t FlowV2(TTree *kine, Double_t psi) {
148 Int_t mult = kine->GetEntries();
152 Double_t Ssin = 0, Scos = 0;
154 kine->Draw("Particles->Phi()", "", "goff");
157 for (Int_t i=0; i<mult; i++) {
159 Ssin += sin( 2 * (phi[i] - psi) );
160 Scos += cos( 2 * (phi[i] - psi) );
163 V2 = sqrt(Ssin*Ssin + Scos*Scos) / mult;
168 ////////////////////////////////////////////////////////////////////////////////