1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // This class analyses TPC cosmics data from the ESD and the ESDfriend
21 // Authors: Jan.Fiete.Grosse-Oetringhaus@cern.ch, Claus.Jorgensen@cern.ch
24 #include "AliROCESDAnalysisSelector.h"
28 #include <AliESDfriend.h>
29 #include <AliTPCclusterMI.h>
30 #include <AliTPCseed.h>
37 #include <TObjArray.h>
38 #include <TTimeStamp.h>
40 #include "TPC/AliTPCClusterHistograms.h"
42 extern TSystem* gSystem;
44 ClassImp(AliROCESDAnalysisSelector)
46 AliROCESDAnalysisSelector::AliROCESDAnalysisSelector() :
52 // Constructor. Initialization of pointers
54 fMinNumberOfRowsIsTrack = 5;
56 fObjectsToSave = new TObjArray();
58 for (Int_t i=0; i<kTPCHists; i++)
59 fClusterHistograms[i] = 0;
62 AliROCESDAnalysisSelector::~AliROCESDAnalysisSelector()
69 void AliROCESDAnalysisSelector::SlaveBegin(TTree* tree)
73 AliSelector::SlaveBegin(tree);
76 void AliROCESDAnalysisSelector::Init(TTree *tree)
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.
84 AliSelector::Init(tree);
86 printf("Init called %p\n", (void*) fESDfriend);
91 tree->SetBranchAddress("ESDfriend", &fESDfriend);
93 tree->SetBranchStatus("*", 0);
94 tree->SetBranchStatus("fTracks.*", 1);
95 tree->SetBranchStatus("fTimeStamp", 1);
96 //tree->SetBranchStatus("fTracks.fCalibContainer", 0);
100 AliDebug(AliLog::kInfo, "INFO: Found ESDfriend branch in chain.");
103 Bool_t AliROCESDAnalysisSelector::Process(Long64_t entry)
106 // Implement your analysis here. Do not forget to call the parent class Process by
107 // if (AliSelector::Process(entry) == kFALSE)
111 if (AliSelector::Process(entry) == kFALSE)
114 // Check prerequisites
117 AliDebug(AliLog::kError, "ESD branch not available");
121 // Check prerequisites
124 AliDebug(AliLog::kError, "ESDfriend branch not available");
128 fESD->SetESDfriend(fESDfriend);
130 Int_t flag = ProcessEvent(entry, kFALSE);
132 ProcessEvent(entry, kTRUE);
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....
148 Int_t AliROCESDAnalysisSelector::ProcessEvent(Long64_t entry, Bool_t detailedHistogram)
151 // Looping over tracks and clusters in event and filling histograms
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
158 // save maximum 50 objects
159 if (detailedHistogram)
160 if (fObjectsToSave->GetEntries() > 50)
163 // for saving single events
164 AliTPCClusterHistograms* clusterHistograms[kTPCSectors];
165 for (Int_t i=0; i<kTPCSectors; i++)
166 clusterHistograms[i] = 0;
168 Bool_t intToReturn = 0;
170 Int_t nTracks = fESD->GetNumberOfTracks();
172 Int_t nSkippedSeeds = 0;
173 Int_t nSkippedTracks = 0;
175 // for "flash" detection
177 Float_t clusterQtotSumVsTime[250];
178 for (Int_t z=0; z<250; z++)
179 clusterQtotSumVsTime[z] = 0;
181 // loop over esd tracks
182 for (Int_t t=0; t<nTracks; t++)
184 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (fESD->GetTrack(t));
187 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
191 AliESDfriendTrack* friendtrack = const_cast<AliESDfriendTrack*> (dynamic_cast<const AliESDfriendTrack*> (esdTrack->GetFriendTrack()));
194 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve friend of track %d.", t));
198 const AliTPCseed* seed = dynamic_cast<const AliTPCseed*> (friendtrack->GetCalibObject(0));
201 AliDebug(AliLog::kDebug, Form("ERROR: Could not retrieve seed of track %d.", t));
206 if (!AcceptTrack(seed, fMinNumberOfRowsIsTrack))
208 AliDebug(AliLog::kDebug, Form("INFO: Rejected track %d.", t));
213 for (Int_t clusterID = 0; clusterID < 160; clusterID++)
215 AliTPCclusterMI* cluster = seed->GetClusterPointer(clusterID);
219 Int_t detector = cluster->GetDetector();
221 if (detector < 0 || detector >= kTPCSectors)
223 AliDebug(AliLog::kDebug, Form("We found a cluster from invalid sector %d", detector));
227 if (!detailedHistogram) {
228 // TODO: find a clever way to handle the time
231 if (fESD->GetTimeStamp()>1160000000)
232 time = fESD->GetTimeStamp();
234 if (!fClusterHistograms[detector])
235 fClusterHistograms[detector] = new AliTPCClusterHistograms(detector,"",time,time+5*60*60);
237 if (!fClusterHistograms[detector+kTPCSectors])
238 fClusterHistograms[detector+kTPCSectors] = new AliTPCClusterHistograms(detector,"",time,time+5*60*60, kTRUE);
240 fClusterHistograms[detector]->FillCluster(cluster, time);
241 fClusterHistograms[detector+kTPCSectors]->FillCluster(cluster, time);
243 Int_t z = Int_t(cluster->GetZ());
246 clusterQtotSumVsTime[z] += cluster->GetQ();
248 } // end of if !detailedHistograms
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));
254 clusterHistograms[detector]->FillCluster(cluster);
258 for (Int_t i=0; i<kTPCHists; i++)
259 if (fClusterHistograms[i])
260 fClusterHistograms[i]->FillTrack(seed);
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) {
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()));
274 for (Int_t i=0; i< kTPCSectors; i++) {
275 if (clusterHistograms[i]) {
276 fObjectsToSave->Add(clusterHistograms[i]);
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);
289 Bool_t AliROCESDAnalysisSelector::AcceptTrack(const AliTPCseed* track, Int_t minRowsIncluded) {
291 // check if the track should be accepted.
293 const Int_t kMinClusters = 5;
294 const Float_t kMinRatio = 0.75;
295 const Float_t kMax1pt = 0.5;
297 Int_t nRowsUsedByTracks = 0;
298 Bool_t rowIncluded[96];
300 Float_t totalQtot = 0;
303 for(Int_t r=0; r<96; r++)
304 rowIncluded[r] = kFALSE;
306 for (Int_t clusterID = 0; clusterID < 160; clusterID++) {
307 AliTPCclusterMI* cluster = track->GetClusterPointer(clusterID);
312 Float_t qTot = cluster->GetQ();
317 if (!rowIncluded[cluster->GetRow()]) {
319 rowIncluded[cluster->GetRow()] = kTRUE;
323 Float_t meanQtot = totalQtot/nClusters;
328 if (nRowsUsedByTracks < minRowsIncluded)
331 // printf(Form(" TRACK: n clusters = %d, n pad rows = %d \n",track->GetNumberOfClusters(), nRowsUsedByTracks));
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;
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;
346 void AliROCESDAnalysisSelector::SlaveTerminate()
352 for (Int_t i=0; i<kTPCHists; i++)
353 if (fClusterHistograms[i])
354 fOutput->Add(fClusterHistograms[i]);
358 void AliROCESDAnalysisSelector::Terminate()
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
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)));
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"));
383 AliDebug(AliLog::kInfo, Form("INFO: Found comment in input list: %s", comment->GetTitle()));
388 TFile* file = TFile::Open(Form("rocESD_%s.root",comment->GetTitle()), "RECREATE");
390 for (Int_t i=0; i<kTPCHists; i++)
391 if (fClusterHistograms[i]) {
392 fClusterHistograms[i]->SaveHistograms();
393 TCanvas* c = fClusterHistograms[i]->DrawHistograms();
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()));
404 gDirectory->mkdir("saved_objects");
405 gDirectory->cd("saved_objects");
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();
413 fObjectsToSave->At(i)->Write();
417 gDirectory->cd("../");