2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
19 gSystem->Load("libANALYSIS");
20 gSystem->Load("libTPCcalib");
24 gSystem->Load("libANALYSIS");
25 gSystem->Load("libTPCcalib");
27 TFile f("CalibObjects.root");
28 AliTPCcalibTrigger *calibTrigger = (AliTPCcalibTrigger *)f->Get("TPCCalib")->FindObject("calibTrigger");
31 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
32 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
33 AliXRDPROOFtoolkit tool;
34 TChain * chainTrack = tool.MakeChain("trigger.txt","Track",0,10200);
36 TChain * chainEvent = tool.MakeChain("trigger.txt","Event",0,10200);
42 #include "Riostream.h"
48 #include "THnSparse.h"
56 #include "TGraphErrors.h"
59 #include "AliTPCclusterMI.h"
60 #include "AliTPCseed.h"
61 #include "AliESDVertex.h"
62 #include "AliESDEvent.h"
63 #include "AliESDfriend.h"
64 #include "AliESDInputHandler.h"
65 #include "AliAnalysisManager.h"
67 #include "AliTracker.h"
69 #include "AliTPCCalROC.h"
73 #include "AliTPCcalibTrigger.h"
75 #include "TTreeStream.h"
76 #include "AliTPCTracklet.h"
77 #include "TTimeStamp.h"
78 #include "AliTPCcalibDB.h"
79 #include "AliTPCcalibLaser.h"
80 #include "AliDCSSensorArray.h"
81 #include "AliDCSSensor.h"
83 ClassImp(AliTPCcalibTrigger)
85 AliTPCcalibTrigger::AliTPCcalibTrigger():
86 AliTPCcalibBase("calibTrigger","calibTrigger"),
92 AliTPCcalibTrigger::AliTPCcalibTrigger(const char * name, const char * title):
93 AliTPCcalibBase(name,title),
102 Long64_t AliTPCcalibTrigger::Merge(TCollection *li) {
106 TIterator* iter = li->MakeIterator();
107 AliTPCcalibTrigger* cal = 0;
109 while ((cal = (AliTPCcalibTrigger*)iter->Next())) {
110 if (!cal->InheritsFrom(AliTPCcalibTrigger::Class())) {
111 Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
114 TMap * addMap=(cal->fHisMap);
115 if(!addMap) return 0;
116 TIterator* iterator = addMap->MakeIterator();
121 while((object=iterator->Next())){
122 THnSparse* his1 = dynamic_cast<THnSparseF*>(cal->fHisMap->GetValue(object->GetName()));
125 THnSparse* his0 = dynamic_cast<THnSparseF*>(fHisMap->GetValue(object->GetName()));
128 his0=MakeHisto(object->GetName());
129 fHisMap->Add(new TObjString(object->GetName()),his0);
139 void AliTPCcalibTrigger::Process(AliESDEvent *event){
144 const TString &trigger = event->GetFiredTriggerClasses();
145 TTreeSRedirector * cstream = GetDebugStreamer();
147 TObjString str(event->GetFiredTriggerClasses());
148 Bool_t hasPIXEL=HasPIXEL(&str);
149 Bool_t hasTRD=HasTRD(&str);
150 Bool_t hasTOF=HasTOF(&str);
151 Bool_t hasACORDE=HasACORDE(&str);
153 if (!GetHisto(trigger.Data())){
154 AddHisto(trigger.Data(),MakeHisto(trigger.Data()));
156 if (!GetHisto("all")){
157 AddHisto("all",MakeHisto("all"));
160 THnSparse *histoAll = GetHisto("all");
161 THnSparse *histo = GetHisto(trigger.Data());
162 Double_t xcont[9]={0,0,0,0,0,0,0,0,0};
164 Int_t ntracks = event->GetNumberOfTracks();
170 AliESDtrack * longest=0;
172 for (Int_t itrack=0; itrack<ntracks; itrack++){
173 AliESDtrack *track=event->GetTrack(itrack);
174 if (!track) continue;
175 if (track->GetTPCNcls()<=nclmax) continue;
176 nclmax = track->GetTPCNcls();
180 histoAll->Fill(xcont);
183 (*cstream) << "Event" <<
189 "acorde="<<hasACORDE<<
190 "ntracks="<<ntracks<<
195 for (Int_t itrack=0; itrack<ntracks; itrack++){
196 AliESDtrack *track=event->GetTrack(itrack);
197 if (!track) continue;
200 track->GetDZ(0.,0.,0.,event->GetMagneticField(),dca);
201 Bool_t status = track->GetPxPyPz(pxyz);
202 Double_t alpha = TMath::ATan2(pxyz[1],pxyz[0]);
203 xcont[1]=track->GetTPCNcls();
207 xcont[5]=track->GetParameter()[3];
208 xcont[6]=track->Pt();
209 xcont[7]=track->GetTPCsignal();
210 histoAll->Fill(xcont);
215 Double_t mpt = track->GetParameter()[4];
216 (*cstream) << "Track" <<
220 "ntracks="<<ntracks<<
224 "acorde="<<hasACORDE<<
238 THnSparse *AliTPCcalibTrigger::MakeHisto(const char* trigger){
240 // Make event/track histograms
243 // ntracks nclMax dcaR dcaZ alpha theta pt dEdx ev
244 Int_t bins[9] = {50, 40, 20, 20, 18, 25, 25, 25, 2 };
245 //Int_t bins[9] = {50* 20* 25* 25* 18* 25* 25* 25 };
246 Double_t xmin[9] = {0., 0, 0, -250, -3.14, -1.5, 0, 0, -1.};
247 Double_t xmax[9] = {50, 160, 150, 250, 3.14, 1.5, 100, 100, 1.};
248 TString axisName[9]={
259 TString axisTitle[9]={
272 THnSparse *sparse = new THnSparseF(Form("his_%s",trigger), Form("his_%s",trigger), 9, bins, xmin, xmax);
273 for (Int_t iaxis=0; iaxis<9; iaxis++){
274 sparse->GetAxis(iaxis)->SetName(axisName[iaxis]);
275 sparse->GetAxis(iaxis)->SetTitle(axisTitle[iaxis]);
280 THnSparse * AliTPCcalibTrigger::GetHisto(const char *trigger) {
282 // return histogram for given class
283 if (!fHisMap) fHisMap=new TMap;
284 return (THnSparse*) fHisMap->GetValue(trigger);
287 void AliTPCcalibTrigger::AddHisto(const char *trigger, THnSparse *his) {
288 if (!GetHisto(trigger)) {
289 TObjString *pstr = new TObjString(trigger);
290 fHisMap->Add(pstr,his);
294 TTree * AliTPCcalibTrigger::MakeTree(const char * fname){
298 TTreeSRedirector * sred = new TTreeSRedirector(fname);
299 TTreeStream &pcstream = (*sred)<<"Trigger";
302 TIterator* iterator = fHisMap->MakeIterator();
305 while((object=iterator->Next())){
306 MakeTree(pcstream, object->GetName());
309 TFile *f = new TFile(fname);
310 TTree *tree = (TTree*)f->Get("Trigger");
315 void AliTPCcalibTrigger::MakeTree(TTreeStream &pcstream, const char *tname){
317 // TTreeSRedirector * sred = new TTreeSRedirector("trigger.root");
318 // TTreeStream &pcstream = (*sred)<<"Trigger";
320 //AliTPCcalibTrigger *calibTrigger = this;
322 THnSparse * his = GetHisto(tname);
325 Int_t *bins = new Int_t[100];
326 Int_t ndim = his->GetNdimensions();
327 Double_t position[10];
329 TObjString str(tname);
330 Bool_t isAll = str.String().Contains("all");
331 Bool_t hasPIXEL=HasPIXEL(&str);
332 Bool_t hasTRD=HasTRD(&str);
333 Bool_t hasTOF=HasTOF(&str);
334 Bool_t hasACORDE=HasACORDE(&str);
335 for (Long64_t i = 0; i < his->GetNbins(); ++i) {
336 value = his->GetBinContent(i, bins);
337 pcstream<<"val="<<value;
338 pcstream<<"tname.="<<&str;
340 pcstream<<"all="<<isAll;
341 pcstream<<"pixel="<<hasPIXEL;
342 pcstream<<"trd="<<hasTRD;
343 pcstream<<"tof="<<hasTOF;
344 pcstream<<"acorde="<<hasACORDE;
346 for (Int_t idim = 0; idim < ndim; idim++) {
347 position[idim] = his->GetAxis(idim)->GetBinCenter(bins[idim]);
348 pcstream<<Form("%s=",his->GetAxis(idim)->GetName())<<position[idim];
355 Bool_t AliTPCcalibTrigger::HasTOF(TObjString *tname){
357 Bool_t result = kFALSE;
358 result|=(tname->String().Contains("0OB")>0);
359 result|=(tname->String().Contains("0OC")>0);
363 Bool_t AliTPCcalibTrigger::HasACORDE(TObjString *tname){
364 Bool_t result = kFALSE;
365 result|=(tname->String().Contains("0ASL")>0);
366 result|=(tname->String().Contains("0AMU")>0);
367 result|=(tname->String().Contains("0ASC")>0);
371 Bool_t AliTPCcalibTrigger::HasPIXEL(TObjString *tname){
372 return (tname->String().Contains("0SCO")>0);
375 Bool_t AliTPCcalibTrigger::HasTRD(TObjString *tname){
376 return (tname->String().Contains("1H")>0);