PWGPP-71 - Adding ncluster statistic to optional streamers. Print additional inform...
[u/mrichter/AliRoot.git] / TPC / FindKrClusters.C
1
2 /****************************************************************************
3  *           Origin: A.Matyja amatyja@cern.ch                               *
4  ****************************************************************************/
5
6 /*
7
8   macro to create array of clusters from TPC digits
9   input files - galice.root 
10                 digits.root - file with digits - usualy use link to galice.root
11                             - in splitted mode - neccesary to create link to proper file
12                             
13    output file - TPC.RecPoints.root
14
15
16 //  Warning - if cluster file AliTPCclusters.root already exist - macro exit and don't produce anything
17                
18  
19 */
20
21
22 #ifndef __CINT__
23 #include <iostream.h>
24 #include "AliRun.h"
25 #include "AliTPCv4.h"
26 #include "AliTPCParam.h"
27 #include "AliTPCclusterKr.h"
28 #include "AliTPCclustererKr.h"
29 #include "TFile.h"
30 #include "TStopwatch.h"
31 #include "TTree.h"
32 #endif
33
34 Int_t FindKrClusters(){
35
36   //
37   //Load DataBase
38   //
39   //char *ocdbpath ="local:///afs/cern.ch/alice/tpctest/OCDB";
40   //char *ocdbpath ="local:///home/matyja/baza/OCDB";
41   char *ocdbpath ="local:///data/baza/OCDB";
42   if (ocdbpath==0){
43     ocdbpath="alien://folder=/alice/data/2007/LHC07w/OCDB/";
44   }
45   printf("OCDB PATH = %s\n",ocdbpath); 
46   AliCDBManager * man = AliCDBManager::Instance();
47   man->SetDefaultStorage(ocdbpath);
48   man->SetRun(0);
49
50   AliRunLoader* rl = AliRunLoader::Open("galice.root");
51   if (rl == 0x0) {
52     cerr<<"Can not open session"<<endl;
53     return 1;
54   }
55   
56   AliTPCLoader *tpcl = (AliTPCLoader*)rl->GetLoader("TPCLoader");
57   if (tpcl == 0x0) {
58     cerr<<"Can not get TPC Loader"<<endl;
59     return 1;
60   }
61
62   if (tpcl->LoadDigits()) {
63     cerr<<"Error occured while loading digits"<<endl;
64     return 1;
65   }
66
67   if (rl->LoadgAlice()) {
68     cerr<<"Error occured while LoadgAlice"<<endl;
69     return 1;
70   }
71   
72   gAlice=rl->GetAliRun();
73   if (!gAlice) {
74     cerr<<"Can't get gAlice !\n";
75     return 1;
76   }
77
78   TDirectory *cwd = gDirectory;
79
80   AliTPCv4 *tpc = (AliTPCv4*)gAlice->GetDetector("TPC");
81   Int_t ver = tpc->IsVersion(); 
82   cerr<<"TPC version "<<ver<<" has been found !\n";
83
84   rl->CdGAFile();
85   
86   AliTPCParam *param=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60");
87   if (!param) {cerr<<"TPC parameters have not been found !\n"; return 4;}
88   
89   AliTPCDigitsArray *digarr=new AliTPCDigitsArray;
90   digarr->Setup(param);
91
92   cerr<<"It has begun"<<endl;  
93   TStopwatch timer;
94   timer.Start();
95   cwd->cd();
96
97   TTree *output_tree;
98   
99   AliTPCclustererKr *clusters = new AliTPCclustererKr();
100   clusters->SetParam(param);
101   clusters->SetOutput(output_tree);
102
103   clusters->SetMinAdc(3);//signal threshold (everything below is treated as 0)
104   clusters->SetMinTimeBins(2);//number of neighbouring timebins
105   clusters->SetMaxPadRangeCm(5.);//distance of the cluster center to the center of a pad (in cm)
106   clusters->SetMaxRowRangeCm(5.);//distance of the cluster center to the center of a padrow (in cm)
107   clusters->SetMaxTimeRange(7.);//distance of the cluster center to the max time bin on a pad (in tackts)
108   //ie. fabs(centerT - time)<7
109
110   clusters->SetIsolCut(3);//set isolation cut threshold
111   clusters->SetValueToSize(3.1);//cut reduce peak at 0
112
113   Int_t nevmax=rl->GetNumberOfEvents();//number of events in run
114   for(Int_t nev=0;nev<nevmax /*&& nev<1*/ ;nev++){
115     rl->GetEvent(nev);
116     
117     TTree* input_tree= tpcl->TreeD();//tree with digits
118     if (input_tree == 0x0){
119       cerr << "Can not get TreeD for event " <<nev<<endl;
120       continue;
121     }
122     digarr->ConnectTree(input_tree);
123     clusters->SetInput(input_tree);
124     clusters->SetDigArr(digarr);
125     cout<<"Processing event "<<nev<<endl;
126     clusters->FinderIO();
127
128   }
129   delete clusters;
130   timer.Stop(); timer.Print();
131   
132   delete rl;//cleans everything
133   
134   return 0;
135 }