]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG0/TPC/AliROCESDAnalysisSelector.cxx
fixed a memory leak, cleaned some code...
[u/mrichter/AliRoot.git] / PWG0 / TPC / AliROCESDAnalysisSelector.cxx
CommitLineData
df71af87 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
1d7991a5 18//
19// This class analyses TPC cosmics data from the ESD and the ESDfriend
df71af87 20//
1d7991a5 21// Authors: Jan.Fiete.Grosse-Oetringhaus@cern.ch, Claus.Jorgensen@cern.ch
df71af87 22//
df71af87 23
24#include "AliROCESDAnalysisSelector.h"
25
26#include <AliLog.h>
27#include <AliESD.h>
28#include <AliESDfriend.h>
899625a7 29#include <AliTPCclusterMI.h>
30#include <AliTPCseed.h>
df71af87 31
df71af87 32#include <TFile.h>
9ecad4f3 33#include <TMath.h>
df71af87 34#include <TTree.h>
35#include <TCanvas.h>
899625a7 36#include <TSystem.h>
9ecad4f3 37#include <TObjArray.h>
899625a7 38#include <TTimeStamp.h>
2d9e89d4 39
40#include "TPC/AliTPCClusterHistograms.h"
df71af87 41
899625a7 42extern TSystem* gSystem;
43
df71af87 44ClassImp(AliROCESDAnalysisSelector)
45
46AliROCESDAnalysisSelector::AliROCESDAnalysisSelector() :
47 AliSelector(),
9ecad4f3 48 fESDfriend(0),
49 fObjectsToSave(0)
df71af87 50{
51 //
52 // Constructor. Initialization of pointers
53 //
9ecad4f3 54 fMinNumberOfRowsIsTrack = 5;
55
56 fObjectsToSave = new TObjArray();
9cc7192c 57
1d7991a5 58 for (Int_t i=0; i<kTPCHists; i++)
9cc7192c 59 fClusterHistograms[i] = 0;
df71af87 60}
61
62AliROCESDAnalysisSelector::~AliROCESDAnalysisSelector()
63{
64 //
65 // Destructor
66 //
67}
68
69void AliROCESDAnalysisSelector::SlaveBegin(TTree* tree)
70{
71 //
72
73 AliSelector::SlaveBegin(tree);
df71af87 74}
75
76void AliROCESDAnalysisSelector::Init(TTree *tree)
77{
78 // The Init() function is called when the selector needs to initialize
79 // a new tree or chain. Typically here the branch addresses of the tree
80 // will be set. It is normaly not necessary to make changes to the
81 // generated code, but the routine can be extended by the user if needed.
82 // Init() will be called many times when running with PROOF.
83
84 AliSelector::Init(tree);
1d7991a5 85
86 printf("Init called %p\n", (void*) fESDfriend);
df71af87 87
88 // Set branch address
1d7991a5 89 if (tree)
90 {
91 tree->SetBranchAddress("ESDfriend", &fESDfriend);
92
6c3d6245 93 tree->SetBranchStatus("*", 0);
94 tree->SetBranchStatus("fTracks.*", 1);
95 tree->SetBranchStatus("fTimeStamp", 1);
1d7991a5 96 //tree->SetBranchStatus("fTracks.fCalibContainer", 0);
6c3d6245 97 }
1d7991a5 98
df71af87 99 if (fESDfriend != 0)
100 AliDebug(AliLog::kInfo, "INFO: Found ESDfriend branch in chain.");
101}
102
103Bool_t AliROCESDAnalysisSelector::Process(Long64_t entry)
104{
105 //
106 // Implement your analysis here. Do not forget to call the parent class Process by
107 // if (AliSelector::Process(entry) == kFALSE)
108 // return kFALSE;
109 //
110
111 if (AliSelector::Process(entry) == kFALSE)
112 return kFALSE;
113
114 // Check prerequisites
115 if (!fESD)
116 {
117 AliDebug(AliLog::kError, "ESD branch not available");
118 return kFALSE;
119 }
120
121 // Check prerequisites
122 if (!fESDfriend)
123 {
124 AliDebug(AliLog::kError, "ESDfriend branch not available");
125 return kFALSE;
126 }
1d7991a5 127
df71af87 128 fESD->SetESDfriend(fESDfriend);
129
899625a7 130 Int_t flag = ProcessEvent(entry, kFALSE);
131 if (flag == 1)
132 ProcessEvent(entry, kTRUE);
9ecad4f3 133
134 // TODO This should not be needed, the TTree::GetEntry() should take care of this, maybe because it has a reference member, to be analyzed
135 // if the ESDfriend is not deleted we get a major memory leak
136 // here the esdfriend seems to be also deleted, very weird behaviour....
137
138 delete fESD;
139 fESD = 0;
140
141 //delete fESDfriend;
142 //fESDfriend = 0;
143
144
145 return kTRUE;
146}
147
899625a7 148Int_t AliROCESDAnalysisSelector::ProcessEvent(Long64_t entry, Bool_t detailedHistogram)
9ecad4f3 149{
150 //
151 // Looping over tracks and clusters in event and filling histograms
152 //
153 // - if detailedHistogram = kTRUE special histograms are saved (in fObjectsToSave)
154 // - the method returns
155 // 1 : if a "flash" is detected something special in this event
156 //
157
899625a7 158 // save maximum 50 objects
9ecad4f3 159 if (detailedHistogram)
899625a7 160 if (fObjectsToSave->GetEntries() > 50)
9ecad4f3 161 return 0;
899625a7 162
9ecad4f3 163 // for saving single events
164 AliTPCClusterHistograms* clusterHistograms[kTPCSectors];
165 for (Int_t i=0; i<kTPCSectors; i++)
166 clusterHistograms[i] = 0;
167
168 Bool_t intToReturn = 0;
169
df71af87 170 Int_t nTracks = fESD->GetNumberOfTracks();
171
9cc7192c 172 Int_t nSkippedSeeds = 0;
9ecad4f3 173 Int_t nSkippedTracks = 0;
174
175 // for "flash" detection
176 Int_t nClusters = 0;
177 Float_t clusterQtotSumVsTime[250];
178 for (Int_t z=0; z<250; z++)
179 clusterQtotSumVsTime[z] = 0;
9cc7192c 180
df71af87 181 // loop over esd tracks
182 for (Int_t t=0; t<nTracks; t++)
183 {
184 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (fESD->GetTrack(t));
185 if (!esdTrack)
186 {
187 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
188 continue;
189 }
190
191 AliESDfriendTrack* friendtrack = const_cast<AliESDfriendTrack*> (dynamic_cast<const AliESDfriendTrack*> (esdTrack->GetFriendTrack()));
192 if (!friendtrack)
193 {
194 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve friend of track %d.", t));
195 continue;
196 }
197
198 const AliTPCseed* seed = dynamic_cast<const AliTPCseed*> (friendtrack->GetCalibObject(0));
199 if (!seed)
200 {
9cc7192c 201 AliDebug(AliLog::kDebug, Form("ERROR: Could not retrieve seed of track %d.", t));
202 nSkippedSeeds++;
df71af87 203 continue;
204 }
205
9ecad4f3 206 if (!AcceptTrack(seed, fMinNumberOfRowsIsTrack))
207 {
208 AliDebug(AliLog::kDebug, Form("INFO: Rejected track %d.", t));
209 nSkippedTracks++;
0b3ccaa2 210 continue;
9ecad4f3 211 }
212
df71af87 213 for (Int_t clusterID = 0; clusterID < 160; clusterID++)
214 {
215 AliTPCclusterMI* cluster = seed->GetClusterPointer(clusterID);
9ecad4f3 216 if (!cluster)
df71af87 217 continue;
df71af87 218
9cc7192c 219 Int_t detector = cluster->GetDetector();
220
899625a7 221 if (detector < 0 || detector >= kTPCSectors)
222 {
223 AliDebug(AliLog::kDebug, Form("We found a cluster from invalid sector %d", detector));
224 continue;
9cc7192c 225 }
483b1eb8 226
9ecad4f3 227 if (!detailedHistogram) {
899625a7 228 // TODO: find a clever way to handle the time
229 Int_t time = 0;
230
231 if (fESD->GetTimeStamp()>1160000000)
232 time = fESD->GetTimeStamp();
233
234 if (!fClusterHistograms[detector])
235 fClusterHistograms[detector] = new AliTPCClusterHistograms(detector,"",time,time+5*60*60);
236
237 if (!fClusterHistograms[detector+kTPCSectors])
238 fClusterHistograms[detector+kTPCSectors] = new AliTPCClusterHistograms(detector,"",time,time+5*60*60, kTRUE);
239
240 fClusterHistograms[detector]->FillCluster(cluster, time);
241 fClusterHistograms[detector+kTPCSectors]->FillCluster(cluster, time);
242
243 Int_t z = Int_t(cluster->GetZ());
244 if (z>=0 && z<250) {
245 nClusters++;
246 clusterQtotSumVsTime[z] += cluster->GetQ();
247 }
9ecad4f3 248 } // end of if !detailedHistograms
249 else {
899625a7 250 // if we need the detailed histograms for this event
251 if (!clusterHistograms[detector])
252 clusterHistograms[detector] = new AliTPCClusterHistograms(detector, Form("flash_entry%d", entry));
253
254 clusterHistograms[detector]->FillCluster(cluster);
9ecad4f3 255 }
df71af87 256 }
9ecad4f3 257
258 for (Int_t i=0; i<kTPCHists; i++)
259 if (fClusterHistograms[i])
899625a7 260 fClusterHistograms[i]->FillTrack(seed);
261
df71af87 262 }
9cc7192c 263
9ecad4f3 264 // check if there's a very large q deposit ("flash")
265 if (!detailedHistogram) {
266 for (Int_t z=0; z<250; z++) {
267 if (clusterQtotSumVsTime[z] > 150000) {
899625a7 268 printf(Form(" \n -> Entry %lld sum of clusters at time %d is %f, ESD timestamp: %s (%d) \n \n", entry, z, clusterQtotSumVsTime[z], TTimeStamp(fESD->GetTimeStamp()).AsString(), fESD->GetTimeStamp()));
269 intToReturn = 1;
9ecad4f3 270 }
271 }
272 }
273 else {
274 for (Int_t i=0; i< kTPCSectors; i++) {
275 if (clusterHistograms[i]) {
899625a7 276 fObjectsToSave->Add(clusterHistograms[i]);
9ecad4f3 277 }
278 }
279 }
1d7991a5 280
9ecad4f3 281// if (nSkippedSeeds > 0)
282// printf("WARNING: The seed was not found for %d out of %d tracks.\n", nSkippedSeeds, nTracks);
283// if (nSkippedTracks > 0)
284// printf("INFO: Rejected %d out of %d tracks.\n", nSkippedTracks, nTracks);
df71af87 285
9ecad4f3 286 return intToReturn;
df71af87 287}
288
9ecad4f3 289Bool_t AliROCESDAnalysisSelector::AcceptTrack(const AliTPCseed* track, Int_t minRowsIncluded) {
0b3ccaa2 290 //
9ecad4f3 291 // check if the track should be accepted.
0b3ccaa2 292 //
9ecad4f3 293 const Int_t kMinClusters = 5;
0b3ccaa2 294 const Float_t kMinRatio = 0.75;
295 const Float_t kMax1pt = 0.5;
296
9ecad4f3 297 Int_t nRowsUsedByTracks = 0;
298 Bool_t rowIncluded[96];
299
300 Float_t totalQtot = 0;
301 Int_t nClusters = 0;
302
303 for(Int_t r=0; r<96; r++)
304 rowIncluded[r] = kFALSE;
305
306 for (Int_t clusterID = 0; clusterID < 160; clusterID++) {
307 AliTPCclusterMI* cluster = track->GetClusterPointer(clusterID);
308
309 if (!cluster)
310 continue;
311
312 Float_t qTot = cluster->GetQ();
313
314 nClusters++;
315 totalQtot += qTot;
316
317 if (!rowIncluded[cluster->GetRow()]) {
318 nRowsUsedByTracks++;
319 rowIncluded[cluster->GetRow()] = kTRUE;
320 }
321 }
322
323 Float_t meanQtot = totalQtot/nClusters;
324
325 if (meanQtot<70)
326 return kFALSE;
0b3ccaa2 327
9ecad4f3 328 if (nRowsUsedByTracks < minRowsIncluded)
329 return kFALSE;
330
331 // printf(Form(" TRACK: n clusters = %d, n pad rows = %d \n",track->GetNumberOfClusters(), nRowsUsedByTracks));
332
0b3ccaa2 333 if (track->GetNumberOfClusters()<kMinClusters) return kFALSE;
334 Float_t ratio = track->GetNumberOfClusters()/(track->GetNFoundable()+1.);
335 if (ratio<kMinRatio) return kFALSE;
336 Float_t mpt = track->Get1Pt();
337 if (TMath::Abs(mpt)>kMax1pt) return kFALSE;
338
339 //if (TMath::Abs(track->GetZ())>240.) return kFALSE;
340 //if (TMath::Abs(track->GetZ())<10.) return kFALSE;
341 //if (TMath::Abs(track->GetTgl())>0.03) return kFALSE;
342
343 return kTRUE;
344}
345
df71af87 346void AliROCESDAnalysisSelector::SlaveTerminate()
347{
348 //
349
1d7991a5 350 if (fOutput)
351 {
352 for (Int_t i=0; i<kTPCHists; i++)
353 if (fClusterHistograms[i])
9cc7192c 354 fOutput->Add(fClusterHistograms[i]);
1d7991a5 355 }
df71af87 356}
357
358void AliROCESDAnalysisSelector::Terminate()
359{
1d7991a5 360 //
361 // read the objects from the output list and write them to a file
362 // the filename is modified by the object comment passed in the tree info or input list
363 //
364
365 if (fOutput)
366 {
367 fOutput->Print();
368
369 for (Int_t i=0; i<kTPCSectors; i++)
370 fClusterHistograms[i] = dynamic_cast<AliTPCClusterHistograms*> (fOutput->FindObject(AliTPCClusterHistograms::FormDetectorName(i, kFALSE)));
371 for (Int_t i=0; i<kTPCSectors; i++)
372 fClusterHistograms[kTPCSectors+i] = dynamic_cast<AliTPCClusterHistograms*> (fOutput->FindObject(AliTPCClusterHistograms::FormDetectorName(i, kTRUE)));
373 }
374
375 TNamed* comment = 0;
376 if (fTree && fTree->GetUserInfo())
377 comment = dynamic_cast<TNamed*>(fTree->GetUserInfo()->FindObject("comment"));
378 if (!comment && fInput)
379 comment = dynamic_cast<TNamed*>(fInput->FindObject("comment"));
6c3d6245 380
381 if (comment)
1d7991a5 382 {
899625a7 383 AliDebug(AliLog::kInfo, Form("INFO: Found comment in input list: %s", comment->GetTitle()));
1d7991a5 384 }
385 else
386 return;
6c3d6245 387
388 TFile* file = TFile::Open(Form("rocESD_%s.root",comment->GetTitle()), "RECREATE");
899625a7 389
1d7991a5 390 for (Int_t i=0; i<kTPCHists; i++)
483b1eb8 391 if (fClusterHistograms[i]) {
9cc7192c 392 fClusterHistograms[i]->SaveHistograms();
483b1eb8 393 TCanvas* c = fClusterHistograms[i]->DrawHistograms();
899625a7 394 TString dir;
395 dir.Form("WWW/%s/%s", comment->GetTitle(), c->GetName());
396 gSystem->mkdir(dir, kTRUE);
397 c->SaveAs(Form("%s/plots_%s_%s.eps",dir.Data(),comment->GetTitle(),c->GetName()));
398 c->SaveAs(Form("%s/plots_%s_%s.gif",dir.Data(),comment->GetTitle(),c->GetName()));
6c3d6245 399
400 c->Close();
401 delete c;
483b1eb8 402 }
9ecad4f3 403
404 gDirectory->mkdir("saved_objects");
405 gDirectory->cd("saved_objects");
406
407 for (Int_t i=0; i<fObjectsToSave->GetSize(); i++) {
408 if (fObjectsToSave->At(i)) {
409 AliTPCClusterHistograms* clusterHistograms = dynamic_cast<AliTPCClusterHistograms*> (fObjectsToSave->At(i));
410 if (clusterHistograms)
411 clusterHistograms->SaveHistograms();
412 else
413 fObjectsToSave->At(i)->Write();
414 }
415 }
416
417 gDirectory->cd("../");
418
419
df71af87 420 file->Close();
899625a7 421}