First version of kdtree (Alexander, Marian)
[u/mrichter/AliRoot.git] / MONITOR / AliMonitorHLT.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 //  This class creates and fills the monitor histograms for the HLT          //
21 //                                                                           //
22 ///////////////////////////////////////////////////////////////////////////////
23
24
25 #include "AliMonitorHLT.h"
26 #include "AliMonitorTrend.h"
27 #include "AliTPCParam.h"
28 #include "AliLog.h"
29 #include <TFolder.h>
30 #include <stdlib.h>
31 #include "AliHLTMemHandler.h"
32 #include "AliHLTSpacePointData.h"
33 #include "AliHLTTrackArray.h"
34 #include "AliHLTTrack.h"
35 #include "AliHLTTransform.h"
36 #include "AliHLTVertex.h"
37
38 //_____________________________________________________________________________
39 AliMonitorHLT::AliMonitorHLT(AliTPCParam* param):
40   AliMonitor(),
41   fParam(param),
42   fClustersCharge(NULL),
43   fNClustersVsRow(NULL),
44   fNClustersVsSector(NULL),
45   fNTracks(NULL),
46   fTrackPt(NULL),
47   fTrackEta(NULL),
48   fTrackPhi(NULL),
49   fTrackNHits(NULL),
50   fTrackDEdxVsP(NULL),
51   fTrackDEdx(NULL),
52   fTrackDz0(NULL),
53   fTrackDr0(NULL),
54   fTrackEtaVsPhi(NULL),
55   fPtEtaVsPhi(NULL),
56   fTrackZvsNHits(NULL),
57   fTrackXYvsNHits(NULL)
58 {
59 // create a HLT monitor object with the given parameters
60
61 }
62
63 //_____________________________________________________________________________
64 void AliMonitorHLT::CreateHistos(TFolder* folder)
65 {
66 // create the HLT monitor histograms
67
68   fFolder = folder->AddFolder("HLT", "HLT");
69
70   fClustersCharge = CreateHisto1("ClustersCharge", 
71                                  "charge distribution of clusters", 
72                                  100, 0, 1000, "charge", "#Delta N/N",
73                                  AliMonitorHisto::kNormEvents);
74
75   Int_t nRows = fParam->GetNRowLow() + fParam->GetNRowUp();
76   fNClustersVsRow = CreateHisto1("NClustersVsRow", 
77                                  "mean number of clusters per pad row", 
78                                  nRows, -0.5, nRows-0.5,
79                                  "pad row", "<N_{clusters}>",
80                                  AliMonitorHisto::kNormEvents);
81
82   Int_t nSector = fParam->GetNInnerSector();
83   fNClustersVsSector = CreateHisto1("NClustersVsSector", 
84                                     "mean number of clusters per sector", 
85                                     nSector, -0.5, nSector-0.5, 
86                                     "sector", "<N_{clusters}>",
87                                     AliMonitorHisto::kNormEvents);
88
89   fNTracks = CreateTrend("NTracks", "number of tracks per event", 
90                          "N_{tracks}");
91
92   fTrackPt = CreateHisto1("TrackPt", "pt distribution of tracks", 
93                           90, 0, 3, "p_{t} [GeV/c]", "#Delta N/N",
94                           AliMonitorHisto::kNormNone);
95
96   fTrackEta = CreateHisto1("TrackEta", "eta distribution of tracks", 
97                            100, -2, 2, "#eta", "#Delta N/N",
98                            AliMonitorHisto::kNormEntries);
99
100   fTrackPhi = CreateHisto1("TrackPhi", "phi distribution of tracks", 
101                            120, 0, 360, "#phi [#circ]", "#Delta N/N",
102                            AliMonitorHisto::kNormEntries);
103
104   fTrackNHits = CreateHisto1("TrackNHits", "Number of hits per track", 
105                            200, 0, 200, "N_{hits}", "#Delta N/N",
106                            AliMonitorHisto::kNormNone);
107
108   fTrackXYvsNHits = CreateHisto2("TrackXYvsNHits", "XY vs Number of hits per track", 
109                                  50, 0, 200,50,0,10,
110                                  "N_{hits}","Radius [cm]","#Delta N/N",
111                                  AliMonitorHisto::kNormNone);
112
113   fTrackZvsNHits = CreateHisto2("TrackZvsNHits", "Z vs Number of hits per track", 
114                                  50, 0, 200,50,-20,20,
115                                  "N_{hits}","Z [cm]","#Delta N/N",
116                                  AliMonitorHisto::kNormNone);
117
118   fTrackDEdxVsP = CreateHisto2("TrackDEdxVsP", "dE/dx of tracks", 
119                                100, 0, 3, 100, 0, 1000, 
120                                "p [GeV/c]", "dE/dx", "#Delta N/N",
121                                AliMonitorHisto::kNormEntries);
122
123   fTrackDEdx = CreateHisto1("TrackDEdx", "dE/dx for tracks with 0.4<p<1.0 GeV/c", 
124                                50, 0, 300, 
125                                "dE/dx", "#Delta N/N",
126                                AliMonitorHisto::kNormEntries);
127
128   fTrackDz0 = CreateHisto1("TrackDz0", "Dz0 of tracks", 
129                            100, -100, 100, "#Delta z0 [cm]", "#Delta N/N",
130                            AliMonitorHisto::kNormEntries);
131
132   fTrackDr0 = CreateHisto1("TrackDr0", "Dr0 of tracks", 
133                            130, 80, 250, "#Delta r0 [cm]", "#Delta N/N",
134                            AliMonitorHisto::kNormEntries);
135
136   fTrackEtaVsPhi = CreateHisto2("TrackEtaVsPhi", "#phi vs #eta", 
137                                20, -1, 1, 25, 0, 360, 
138                                "#eta", "#phi", "#Delta N/N",
139                                AliMonitorHisto::kNormNone);
140
141   fPtEtaVsPhi = CreateHisto2("PtEtaVsPhi", "#phi vs #eta", 
142                                20, -1, 1, 25, 0, 360, 
143                                "#eta", "#phi", "#Delta N/N",
144                                AliMonitorHisto::kNormNone);
145
146 }
147
148
149 //_____________________________________________________________________________
150 void AliMonitorHLT::FillHistos(AliRunLoader* /*runLoader*/, 
151                                AliRawReader* /*rawReader*/, AliESDEvent* /*esd*/)
152 {
153 // fill the HLT monitor histogrms
154
155   AliHLTMemHandler clusterHandler[36];
156   AliHLTSpacePointData *clusters[36];
157   for (Int_t iSector = 0; iSector < fParam->GetNInnerSector(); iSector++) {
158     char fileName[256];
159     sprintf(fileName, "hlt/points_%d_-1.raw", iSector);
160     if (!clusterHandler[iSector].SetBinaryInput(fileName)) {
161       AliWarning(Form("could not open file %s", fileName));
162       continue;
163     }
164     clusters[iSector] = (AliHLTSpacePointData*) clusterHandler[iSector].Allocate();
165     UInt_t nClusters = 0;
166     clusterHandler[iSector].Binary2Memory(nClusters, clusters[iSector]);
167
168     for (UInt_t iCluster = 0; iCluster < nClusters; iCluster++) {
169       AliHLTSpacePointData& cluster = clusters[iSector][iCluster];
170       fClustersCharge->Fill(cluster.fCharge);
171       fNClustersVsRow->Fill(cluster.fPadRow);
172       fNClustersVsSector->Fill(iSector);
173     }
174
175     clusterHandler[iSector].CloseBinaryInput();
176   }
177
178   fNClustersVsSector->ScaleErrorBy(10.);
179
180   AliHLTMemHandler memHandler;
181   if (!memHandler.SetBinaryInput("hlt/tracks.raw")) {
182     AliWarning("could not open file hlt/tracks.raw");
183     return;
184   }
185   AliHLTTrackArray* tracks = new AliHLTTrackArray;
186   memHandler.Binary2TrackArray(tracks);
187   Double_t xc,yc,zc;
188   AliHLTVertex vertex;
189
190   fNTracks->Fill(tracks->GetNTracks());
191   for (Int_t iTrack = 0; iTrack < tracks->GetNTracks(); iTrack++) {
192     AliHLTTrack* track = tracks->GetCheckedTrack(iTrack);
193     if(!track) continue;
194     track->CalculateHelix();
195     track->GetClosestPoint(&vertex,xc,yc,zc);
196     if(TMath::Abs(zc)>10.) continue;
197     fTrackPt->Fill(track->GetPt());
198     fTrackEta->Fill(track->GetPseudoRapidity());
199     fTrackPhi->Fill(track->GetPsi() * TMath::RadToDeg());
200     if(track->GetPt()>3.) {
201       fTrackEtaVsPhi->Fill(track->GetPseudoRapidity(),track->GetPsi() * TMath::RadToDeg());
202       fPtEtaVsPhi->Fill(track->GetPseudoRapidity(),track->GetPsi() * TMath::RadToDeg(),track->GetPt());
203     }
204     fTrackDz0->Fill(track->GetZ0());
205     fTrackDr0->Fill(track->GetR0());
206     fTrackNHits->Fill(track->GetNHits());
207     fTrackXYvsNHits->Fill(track->GetNHits(),TMath::Sqrt(xc*xc+yc*yc));
208     fTrackZvsNHits->Fill(track->GetNHits(),zc);
209
210     // Track dEdx
211     Int_t nc=track->GetNHits();
212     UInt_t *hits = track->GetHitNumbers();
213     Float_t sampleDEdx[159];
214     for (Int_t iHit = 0; iHit < nc; iHit++) {
215       UInt_t hitID = hits[iHit];
216       Int_t iSector = (hitID>>25) & 0x7f;
217       UInt_t position = hitID&0x3fffff;
218       UChar_t padrow = clusters[iSector][position].fPadRow;
219       Float_t pWidth = AliHLTTransform::GetPadPitchWidthLow();
220       if (padrow>63)
221         pWidth = AliHLTTransform::GetPadPitchWidthUp(); 
222       Float_t corr=1.; if (padrow>63) corr=0.67;
223       sampleDEdx[iHit] = clusters[iSector][position].fCharge/pWidth*corr;
224       Double_t crossingangle = track->GetCrossingAngle(padrow,iSector);
225       Double_t s = sin(crossingangle);
226       Double_t t = track->GetTgl();
227       sampleDEdx[iHit] *= TMath::Sqrt((1-s*s)/(1+t*t));
228     }
229
230     /* Cook dEdx */
231     Int_t i;
232     Int_t swap;//stupid sorting
233     do {
234       swap=0;
235       for (i=0; i<nc-1; i++) {
236         if (sampleDEdx[i]<=sampleDEdx[i+1]) continue;
237         Float_t tmp=sampleDEdx[i];
238         sampleDEdx[i]=sampleDEdx[i+1]; sampleDEdx[i+1]=tmp;
239         swap++;
240       }
241     } while (swap);
242
243     Double_t low=0.05; Double_t up=0.7;
244     Int_t nl=Int_t(low*nc), nu=Int_t(up*nc);
245     Float_t trackDEdx=0;
246     for (i=nl; i<=nu; i++) trackDEdx += sampleDEdx[i];
247     trackDEdx /= (nu-nl+1);
248
249     fTrackDEdxVsP->Fill(track->GetP(),trackDEdx);
250     if(track->GetP()>0.4 && track->GetP()<1.0)
251       fTrackDEdx->Fill(trackDEdx);
252   }
253
254   delete tracks;
255   memHandler.CloseBinaryInput();
256 }