Retrofeed from 4-01-Release
[u/mrichter/AliRoot.git] / RALICE / wa98 / Wa98Convert.cxx
CommitLineData
5c9fc915 1// $Id$
2
3///////////////////////////////////////////////////////////////////////////
4// Class Wa98Convert
5// Conversion of Wa98 Hbook ntuple data into Wa98Event physics event structures.
6//
7// Usage example :
8// ---------------
9//
10// gSystem->Load("libHbook");
11// gSystem->Load("ralice");
12// gSystem->Load("rwa98");
13//
14// // Input file with ANALYZ produced Hbook ntuple data
15// THbookFile* f=new THbookFile("pb611258.cwn");
16// TTree* h999=(TTree*)f->Get(999);
17//
18// // Output file for the event structures
19// TFile* ofile=new TFile("run11258.root","RECREATE","WA98 data in RALICE event structure");
20// TTree* otree=new TTree("T","Data of the 1996 Pb+Pb run");
21//
22// Int_t nentries=h999->GetEntries();
23// cout << " Number of entries available : " << nentries << endl;
24// cout << endl;
25//
26// // Limit the number of entries for testing
27// nentries=300;
28//
29// // Print frequency to produce a short summary print every printfreq events
30// Int_t printfreq=10;
31//
32// Convert q(h999);
33// q.Loop(otree,nentries,printfreq);
34//
35// otree->Print();
36//
37// // Close output file
38// ofile->Write();
39// ofile->Close();
40//
41//--- Author: Nick van Eijndhoven 06-jul-2004 Utrecht University
42//- Modified: NvE $Date$ Utrecht University
43///////////////////////////////////////////////////////////////////////////
44
45#include "Wa98Convert.h"
46#include "Riostream.h"
47
48ClassImp(Wa98Convert) // Class implementation to enable ROOT I/O
49
50Wa98Convert::Wa98Convert(TTree* tree)
51{
52// Default constructor.
53// Initialise the input tree (or chain) to be converted.
54// By default tree=0;
55
56 fChain=tree;
57 if (!fChain) return;
58
59 // Link the variables to the branches of the input tree/chain
60 fChain->SetBranchAddress("Jrun",&Jrun);
61 fChain->SetBranchAddress("Jevt",&Jevt);
62 fChain->SetBranchAddress("Jdate",&Jdate);
63 fChain->SetBranchAddress("Jtime",&Jtime);
64 fChain->SetBranchAddress("Jevid",&Jevid);
65 fChain->SetBranchAddress("Jwscal",&Jwscal);
66 fChain->SetBranchAddress("Itword",&Itword);
67 fChain->SetBranchAddress("Zdc",&Zdc);
68 fChain->SetBranchAddress("Emir",&Emir);
69 fChain->SetBranchAddress("Emire",&Emire);
70 fChain->SetBranchAddress("Emirh",&Emirh);
71 fChain->SetBranchAddress("Etm",&Etm);
72 fChain->SetBranchAddress("Etme",&Etme);
73 fChain->SetBranchAddress("Etmh",&Etmh);
74 fChain->SetBranchAddress("Nmod",&Nmod);
75 fChain->SetBranchAddress("Irowl",Irowl);
76 fChain->SetBranchAddress("Icoll",Icoll);
77 fChain->SetBranchAddress("Adcl",Adcl);
78 fChain->SetBranchAddress("Nclu",&Nclu);
79 fChain->SetBranchAddress("Irowc",Irowc);
80 fChain->SetBranchAddress("Icolc",Icolc);
81 fChain->SetBranchAddress("Adcc",Adcc);
82 fChain->SetBranchAddress("Ncluv",&Ncluv);
83 fChain->SetBranchAddress("Iadccv",Iadccv);
84 fChain->SetBranchAddress("Thetacv",Thetacv);
85 fChain->SetBranchAddress("Phicv",Phicv);
86}
87///////////////////////////////////////////////////////////////////////////
88Wa98Convert::~Wa98Convert()
89{
90// Default destructor.
91}
92///////////////////////////////////////////////////////////////////////////
93void Wa98Convert::Loop(TTree* otree,Int_t nentries,Int_t printfreq)
94{
95// Loop over the specified number of entries and convert the
96// ntuple data into the Wa98Event structure.
97// The output will be written on the output tree specified as "otree".
98// If otree=0, a default standard output tree will be created.
99// If nentries<0 (default) all the entries of the input chain
100// will be processed.
101// Every "printfreq" events a short event summary will be printed.
102// The default value is printfreq=1.
103
104 if (fChain==0) return;
105
106 if (nentries<0) nentries=fChain->GetEntriesFast();
107
108 if (!otree) otree=new TTree("T","Data of the 1996 Pb+Pb run");
109
110 Double_t pi=acos(-1.);
111
112 Double_t me=0.51099890221e-3;
113 Double_t mpi=0.13956995;
114 Double_t mkc=0.493677;
115 Double_t mk0=0.497672;
116 Double_t mp=0.93827231;
117 Double_t mlam=1.115683;
118
119 Wa98Event* evt=new Wa98Event();
120
121 // Branch in the tree for the event structure
122 Int_t split=1;
123 Int_t bsize=32000;
124 otree->Branch("Wa98Event","Wa98Event",&evt,bsize,split);
125
126 // The LEDA specific output data
127 AliCalorimeter* ledaup=new AliCalorimeter(44,144);
128 AliCalorimeter* ledalw=new AliCalorimeter(40,144);
129
130 ledaup->SetName("LedaUp");
131 ledalw->SetName("LedaDown");
132
133 evt->InitLeda(ledaup);
134 evt->InitLeda(ledalw);
135
136 TDatime datim;
137 Float_t pos[3],costh;
138 AliSignal s;
139 s.SetName("CPV signal ADC");
140
141 for (Int_t jentry=0; jentry<nentries; jentry++)
142 {
143 fChain->GetEntry(jentry);
144
145 // Reset the complete Event structure
146 evt->Reset();
147
148 evt->SetRunNumber(Jrun);
149 evt->SetEventNumber(Jevt);
150 datim.Set(Jdate,Jtime);
151 evt->SetDayTime(datim);
152 evt->SetProjectile(207,82,158);
153 evt->SetTarget(207,82,0);
154 evt->SetWeight(Jwscal);
155 evt->SetTrig(Itword);
156 evt->SetZdc(Zdc*1000.);
157 evt->SetMiracE(1000.*Emir,Emire,Emirh);
158 evt->SetMiracEt(Etm,Etme,Etmh);
159
160 ledaup->Reset();
161 ledalw->Reset();
162 // Fill calorimeter with module data
163 for (Int_t i=0; i<Nmod; i++)
164 {
165 if (Adcl[i] > 3) // Adc cut of 3 to remove noise
166 {
167 if (Irowl[i] > 0) ledaup->SetSignal(Irowl[i],Icoll[i],Adcl[i]);
168 if (Irowl[i] < 0) ledalw->SetSignal(-Irowl[i],Icoll[i],Adcl[i]);
169 }
170 }
171
172 // Store associated CPV signals
173 for (Int_t j=0; j<Ncluv; j++)
174 {
175 s.Reset();
176 s.SetSignal(Iadccv[j]);
177 pos[1]=Thetacv[j]*pi/180.;
178 pos[2]=Phicv[j]*pi/180.;
179 costh=cos(pos[1]);
180 pos[0]=0;
181 if (costh) pos[0]=2103./costh;
182 s.SetPosition(pos,"sph");
183 pos[0]=0.4;
184 pos[1]=2.2;
185 pos[2]=0;
186 s.SetPositionErrors(pos,"car");
187 if (Phicv[j]>=0. && Phicv[j]<=180.)
188 {
189 ledaup->AddVetoSignal(s);
190 }
191 else
192 {
193 ledalw->AddVetoSignal(s);
194 }
195 }
196
197 evt->AddDevice(ledaup);
198 evt->AddDevice(ledalw);
199
200 if (!(jentry%printfreq))
201 {
202// cout << " Entry : " << jentry << " Run : " << Jrun << " Event : " << Jevt
203// << " Itword : " << Itword << " Etm : " << Etm << endl;
204// cout << " Jdate : " << Jdate << " Jtime : " << Jtime << endl;
205 cout << " Itword : " << Itword << " Nmod : " << Nmod << " Ncluv : " << Ncluv << endl;
206 evt->HeaderData();
207 }
208
209 // Write the complete structure to the output Tree
210 otree->Fill();
211 }
212
213 if (evt) delete evt;
214 if (ledaup) delete ledaup;
215 if (ledalw) delete ledalw;
216}
217///////////////////////////////////////////////////////////////////////////