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