Macro to read PMD.RecPoints.root
[u/mrichter/AliRoot.git] / PMD / AliPMDRecpointRead.C
CommitLineData
2d3c3da3 1// ----------------------------------------------------//
2// //
3// This macro reads the PMD clusters which //
4// are stored in the file "PMD.RecPoints.root" //
5// //
6// ----------------------------------------------------//
7
8#include <Riostream.h>
9#include "TBranch.h"
10#include "TStopwatch.h"
11
12extern AliRun *gAlice;
13
14Int_t AliPMDRecpointRead(Int_t nevent = 1)
15{
16 if (gAlice)
17 {
18 delete gAlice->GetRunLoader();
19 delete gAlice;//if everything was OK here it is already NULL
20 gAlice = 0x0;
21 }
22 AliRunLoader *fRunLoader = AliRunLoader::Open("galice.root","Event","update");
23 if (!fRunLoader)
24 {
25 cerr<<"Can't load RunLoader"<<endl;
26 return 1;
27 }
28 AliLoader *pmdloader = fRunLoader->GetLoader("PMDLoader");
29 Int_t nevent = fRunLoader->GetNumberOfEvents();
30 cout << " * *********** nevent = " << nevent << endl;
31
32 if (pmdloader == 0x0)
33 {
34 cerr<<" ===> Can not find PMD or PMDLoader <===\n";
35 delete fRunLoader;
36 return 2;
37 }
38
39 pmdloader->LoadRecPoints("READ");
40 TClonesArray *fRecpoints;
41 AliPMDUtility *cc = new AliPMDUtility();
42 TH2F *h2 = new TH2F("h2"," ",100,-100.,100.,100,-100.,100.);
43
44 FILE *fpw = fopen("junk_rec.dat","w");
45
46 for (Int_t ievt = 0; ievt < nevent; ievt++)
47 {
48 fRunLoader->GetEvent(ievt);
49 TTree *treeR = pmdloader->TreeR();
50 if (treeR == 0x0)
51 {
52 cout << " Can not get TreeR" << endl;
53 return 3;
54 }
55
56 AliPMDrecpoint1 *pmdrecpoint;
57 TBranch *branch = treeR->GetBranch("PMDRecpoint");
58 branch->SetAddress(&fRecpoints);
59 /**********************************************************************
60 * det : Detector, 0: PRE & 1:CPV *
61 * smn : Serial Module Number from which Super Module Number *
62 * and Unit Module Numbers are extracted *
63 * xpos : x-position of the cluster *
64 * ypos : y-position of the cluster *
65 * THESE xpos & ypos are not the true xpos and ypos *
66 * for some of the unit modules. They are rotated. *
67 * adc : ADC contained in the cluster *
68 * ncell : Number of cells contained in the cluster *
69 * rad : radius of the cluster (1d fit) *
70 * ism : Supermodule number extracted from smn *
71 * ium : Unit module number extracted from smn *
72 * xpad : TRUE x-position of the cluster *
73 * ypad : TRUE y-position of the cluster *
74 **********************************************************************/
75 Int_t ism, ium;
76 Int_t det,smn;
77 Float_t xpos,ypos, xpad, ypad;
78 Float_t adc, ncell, rad;
79 Float_t xx, yy;
80 Int_t nmodules = treeR->GetEntries();
81 cout << " nmodules = " << nmodules << endl;
82 for (Int_t imodule = 0; imodule < nmodules; imodule++)
83 {
84 treeR->GetEntry(imodule);
85 Int_t nentries = fRecpoints->GetLast();
86 for(Int_t ient = 0; ient < nentries+1; ient++)
87 {
88 pmdrecpoint = (AliPMDrecpoint1*)fRecpoints->UncheckedAt(ient);
89 det = (Int_t) pmdrecpoint->GetDetector();
90 smn = (Int_t) pmdrecpoint->GetSMNumber();
91 xpos = pmdrecpoint->GetClusX();
92 ypos = pmdrecpoint->GetClusY();
93 adc = pmdrecpoint->GetClusADC();
94 ncell = pmdrecpoint->GetClusCells();
95 rad = pmdrecpoint->GetClusRadius();
96 //
97 // Now change the xpos and ypos to its original values
98 // for the unit modules which are earlier changed.
99 // xpad and ypad are the real positions.
100 //
101 if(det == 0 || det == 1)
102 {
103 if(smn < 12)
104 {
105 ism = smn/6;
106 ium = smn - ism*6;
107 xpad = ypos;
108 ypad = xpos;
109 }
110 else if( smn >= 12 && smn < 24)
111 {
112 ism = smn/6;
113 ium = smn - ism*6;
114 xpad = xpos;
115 ypad = ypos;
116 }
117 }
118 //
119 // User has to plug in his analysis code here
120 //
121
122 fprintf(fpw,"%d %d %d %d %f %f %f %f %f\n",
123 det,smn,ism,ium,xpad,ypad,adc,ncell,rad);
124 //
125 // Plot the cluster centroid to see the PMD geometry
126 // using the PMD Utility class
127 //
128 if (det == 1)
129 {
130 // Draw only for PRE plane
131 cc->RectGeomCellPos(ism,ium,xpad,ypad,xx,yy);
132 h2->Fill(xx,yy);
133 }
134
135 //
136 // End of the User code
137 //
138 }
139 }
140
141 }
142
143 h2->Draw();
144
145 fclose(fpw);
146 return 0;
147}
148