Introduction of HLT dE/dx (C.Cheshkov)
[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 #ifdef ALI_HLT
27 #include <stdlib.h>
28 #include "AliL3MemHandler.h"
29 #include "AliL3SpacePointData.h"
30 #include "AliL3TrackArray.h"
31 #include "AliL3Track.h"
32 #include "AliL3Transform.h"
33 #endif
34
35 //_____________________________________________________________________________
36 AliMonitorHLT::AliMonitorHLT(AliTPCParam* param)
37 {
38 // create a HLT monitor object with the given parameters
39
40   fParam = param;
41 }
42
43 //_____________________________________________________________________________
44 AliMonitorHLT::~AliMonitorHLT()
45 {
46 }
47
48
49 //_____________________________________________________________________________
50 void AliMonitorHLT::CreateHistos(TFolder* folder)
51 {
52 // create the HLT monitor histograms
53
54   fFolder = folder->AddFolder("HLT", "HLT");
55
56   fClustersCharge = CreateHisto1("ClustersCharge", 
57                                  "charge distribution of clusters", 
58                                  100, 0, 1000, "charge", "#Delta N/N",
59                                  AliMonitorHisto::kNormEvents);
60
61   Int_t nRows = fParam->GetNRowLow() + fParam->GetNRowUp();
62   fNClustersVsRow = CreateHisto1("NClustersVsRow", 
63                                  "mean number of clusters per pad row", 
64                                  nRows, -0.5, nRows-0.5,
65                                  "pad row", "<N_{clusters}>",
66                                  AliMonitorHisto::kNormEvents);
67
68   Int_t nSector = fParam->GetNInnerSector();
69   fNClustersVsSector = CreateHisto1("NClustersVsSector", 
70                                     "mean number of clusters per sector", 
71                                     nSector, -0.5, nSector-0.5, 
72                                     "sector", "<N_{clusters}>",
73                                     AliMonitorHisto::kNormEvents);
74
75   fNTracks = CreateTrend("NTracks", "number of tracks per event", 
76                          "N_{tracks}");
77
78   fTrackPt = CreateHisto1("TrackPt", "pt distribution of tracks", 
79                           90, 0, 3, "p_{t} [GeV/c]", "#Delta N/N",
80                           AliMonitorHisto::kNormEntries);
81
82   fTrackEta = CreateHisto1("TrackEta", "eta distribution of tracks", 
83                            100, -2, 2, "#eta", "#Delta N/N",
84                            AliMonitorHisto::kNormEntries);
85
86   fTrackPhi = CreateHisto1("TrackPhi", "phi distribution of tracks", 
87                            120, 0, 360, "#phi [#circ]", "#Delta N/N",
88                            AliMonitorHisto::kNormEntries);
89
90   fTrackNHits = CreateHisto1("TrackNHits", "Number of hits per track", 
91                            200, 0, 200, "N_{hits}", "#Delta N/N",
92                            AliMonitorHisto::kNormEntries);
93
94   fTrackDEdxVsP = CreateHisto2("TrackDEdxVsP", "dE/dx of tracks", 
95                                100, 0, 3, 100, 0, 1000, 
96                                "p [GeV/c]", "dE/dx", "#Delta N/N",
97                                AliMonitorHisto::kNormEntries);
98
99   fTrackDz0 = CreateHisto1("TrackDz0", "Dz0 of tracks", 
100                            100, -100, 100, "#Delta z0 [cm]", "#Delta N/N",
101                            AliMonitorHisto::kNormEntries);
102
103   fTrackDr0 = CreateHisto1("TrackDr0", "Dr0 of tracks", 
104                            130, 80, 250, "#Delta r0 [cm]", "#Delta N/N",
105                            AliMonitorHisto::kNormEntries);
106
107   fTrackAngle = CreateHisto1("TrackAngle", "azimutal distribution of tracks", 
108                            100, -4, 4, " ", "#Delta N/N",
109                            AliMonitorHisto::kNormEntries);
110
111 }
112
113
114 #include <TCanvas.h>
115 //_____________________________________________________________________________
116 void AliMonitorHLT::FillHistos(AliRunLoader* /*runLoader*/, 
117                                AliRawReader* /*rawReader*/)
118 {
119 // fill the HLT monitor histogrms
120
121 #ifndef ALI_HLT
122   Warning("FillHistos", "the code was compiled without HLT support");
123
124 #else
125   AliL3MemHandler clusterHandler[36];
126   AliL3SpacePointData *clusters[36];
127   for (Int_t iSector = 0; iSector < fParam->GetNInnerSector(); iSector++) {
128     char fileName[256];
129     sprintf(fileName, "hlt/points_%d_-1.raw", iSector);
130     if (!clusterHandler[iSector].SetBinaryInput(fileName)) {
131       Warning("FillHistos", "could not open file %s", fileName);
132       continue;
133     }
134     clusters[iSector] = (AliL3SpacePointData*) clusterHandler[iSector].Allocate();
135     UInt_t nClusters = 0;
136     clusterHandler[iSector].Binary2Memory(nClusters, clusters[iSector]);
137
138     for (UInt_t iCluster = 0; iCluster < nClusters; iCluster++) {
139       AliL3SpacePointData& cluster = clusters[iSector][iCluster];
140       fClustersCharge->Fill(cluster.fCharge);
141       fNClustersVsRow->Fill(cluster.fPadRow);
142       fNClustersVsSector->Fill(iSector);
143     }
144
145     clusterHandler[iSector].CloseBinaryInput();
146   }
147
148   fNClustersVsSector->ScaleErrorBy(10.);
149
150   AliL3MemHandler memHandler;
151   if (!memHandler.SetBinaryInput("hlt/tracks.raw")) {
152     Warning("FillHistos", "could not open file hlt/tracks.raw");
153     return;
154   }
155   AliL3TrackArray* tracks = new AliL3TrackArray;
156   memHandler.Binary2TrackArray(tracks);
157
158   fNTracks->Fill(tracks->GetNTracks());
159   for (Int_t iTrack = 0; iTrack < tracks->GetNTracks(); iTrack++) {
160     AliL3Track* track = tracks->GetCheckedTrack(iTrack);
161     if(!track) continue;
162     track->CalculateHelix();
163     fTrackPt->Fill(track->GetPt());
164     fTrackEta->Fill(track->GetPseudoRapidity());
165     fTrackPhi->Fill(track->GetPsi() * TMath::RadToDeg());
166     fTrackDz0->Fill(track->GetZ0());
167     fTrackDr0->Fill(track->GetR0());
168     fTrackNHits->Fill(track->GetNHits());
169
170     // Track dEdx
171     Int_t nc=track->GetNHits();
172     UInt_t *hits = track->GetHitNumbers();
173     Float_t sampleDEdx[159];
174     for (Int_t iHit = 0; iHit < nc; iHit++) {
175       UInt_t hitID = hits[iHit];
176       Int_t iSector = (hitID>>25) & 0x7f;
177       UInt_t position = hitID&0x3fffff;
178       UChar_t padrow = clusters[iSector][position].fPadRow;
179       Float_t pWidth = AliL3Transform::GetPadPitchWidthLow();
180       if (padrow>63)
181         pWidth = AliL3Transform::GetPadPitchWidthUp(); 
182       Float_t corr=1.; if (padrow>63) corr=0.67;
183       sampleDEdx[iHit] = clusters[iSector][position].fCharge/pWidth*corr;
184       Double_t crossingangle = track->GetCrossingAngle(padrow,iSector);
185       fTrackAngle->Fill(crossingangle);
186       Double_t s = sin(crossingangle);
187       Double_t t = track->GetTgl();
188       sampleDEdx[iHit] *= TMath::Sqrt((1-s*s)/(1+t*t));
189     }
190
191     /* Cook dEdx */
192     Int_t i;
193     Int_t swap;//stupid sorting
194     do {
195       swap=0;
196       for (i=0; i<nc-1; i++) {
197         if (sampleDEdx[i]<=sampleDEdx[i+1]) continue;
198         Float_t tmp=sampleDEdx[i];
199         sampleDEdx[i]=sampleDEdx[i+1]; sampleDEdx[i+1]=tmp;
200         swap++;
201       }
202     } while (swap);
203
204     Double_t low=0.05; Double_t up=0.7;
205     Int_t nl=Int_t(low*nc), nu=Int_t(up*nc);
206     Float_t trackDEdx=0;
207     for (i=nl; i<=nu; i++) trackDEdx += sampleDEdx[i];
208     trackDEdx /= (nu-nl+1);
209
210     fTrackDEdxVsP->Fill(track->GetP(),trackDEdx);
211   }
212
213   delete tracks;
214   memHandler.CloseBinaryInput();
215 #endif
216 }