1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // This class creates and fills the monitor histograms for the TPC //
22 ///////////////////////////////////////////////////////////////////////////////
30 #include "AliMonitorDataTPC.h"
31 #include "AliMonitorHisto.h"
32 #include "AliMonitorTPC.h"
33 #include "AliMonitorTrend.h"
34 #include "AliRawReader.h"
35 #include "AliRunLoader.h"
36 #include "AliTPCClustersRow.h"
37 #include "AliTPCParam.h"
38 #include "AliTPCRawStream.h"
39 #include "AliTPCclusterMI.h"
41 ClassImp(AliMonitorTPC)
44 //_____________________________________________________________________________
45 AliMonitorTPC::AliMonitorTPC(AliTPCParam* param)
47 // create a TPC monitor object with the given parameters
50 fData = new AliMonitorDataTPC(10000);
53 //_____________________________________________________________________________
54 AliMonitorTPC::AliMonitorTPC(const AliMonitorTPC& monitor) :
57 AliFatal("copy constructor not implemented");
60 //_____________________________________________________________________________
61 AliMonitorTPC& AliMonitorTPC::operator = (const AliMonitorTPC& /*monitor*/)
63 AliFatal("assignment operator not implemented");
67 //_____________________________________________________________________________
68 AliMonitorTPC::~AliMonitorTPC()
74 //_____________________________________________________________________________
75 void AliMonitorTPC::CreateHistos(TFolder* folder)
77 // create the TPC monitor histograms
79 fFolder = folder->AddFolder("TPC", "TPC");
81 fPadsCharge = CreateHisto1("PadsCharge", "charge distribution of pads",
82 100, 0, 200, "charge", "#Delta N/N",
83 AliMonitorHisto::kNormEvents);
85 fClustersCharge = CreateHisto1("ClustersCharge",
86 "charge distribution of clusters",
87 100, 0, 1000, "charge", "#Delta N/N",
88 AliMonitorHisto::kNormEvents);
90 Int_t nRows = fParam->GetNRowLow() + fParam->GetNRowUp();
91 fNClustersVsRow = CreateHisto1("NClustersVsRow",
92 "mean number of clusters per pad row",
93 nRows, -0.5, nRows-0.5,
94 "pad row", "<N_{clusters}>",
95 AliMonitorHisto::kNormEvents);
97 Int_t nSector = fParam->GetNInnerSector();
98 fNClustersVsSector = CreateHisto1("NClustersVsSector",
99 "mean number of clusters per sector",
100 nSector, -0.5, nSector-0.5,
101 "sector", "<N_{clusters}>",
102 AliMonitorHisto::kNormEvents);
104 fNTracks = CreateTrend("NTracks", "number of tracks per event",
107 fTrackPt = CreateHisto1("TrackPt", "pt distribution of tracks",
108 90, 0, 3, "p_{t} [GeV/c]", "#Delta N/N",
109 AliMonitorHisto::kNormNone);
111 fTrackEta = CreateHisto1("TrackEta", "eta distribution of tracks",
112 100, -2, 2, "#eta", "#Delta N/N",
113 AliMonitorHisto::kNormEntries);
115 fTrackPhi = CreateHisto1("TrackPhi", "phi distribution of tracks",
116 120, -180, 180, "#phi [#circ]", "#Delta N/N",
117 AliMonitorHisto::kNormEntries);
118 fTrackPhi->SetDescription("The phi distribution should be flat on average.\nIf it is not flat check for dead TPC sectors.");
120 fTrackNCl = CreateHisto1("TrackNCl", "Number of clusters per track",
121 200, 0, 200, "N_{clusters}", "#Delta N/N",
122 AliMonitorHisto::kNormNone);
124 fTrackDEdxVsP = CreateHisto2("TrackDEdxVsP", "dE/dx of tracks",
125 100, 0, 3, 100, 0, 200,
126 "p [GeV/c]", "dE/dx", "#Delta N/N",
127 AliMonitorHisto::kNormEntries);
129 fTrackDEdx = CreateHisto1("TrackDEdx", "dE/dx of tracks with 0.4<p<1.0 GeV/c",
131 "dE/dx", "#Delta N/N",
132 AliMonitorHisto::kNormEntries);
134 fTrackEtaVsPhi = CreateHisto2("TrackEtaVsPhi", "#phi vs #eta",
135 20, -1, 1, 25, 0, 360,
136 "#eta", "#phi", "#Delta N/N",
137 AliMonitorHisto::kNormNone);
139 fPtEtaVsPhi = CreateHisto2("PtEtaVsPhi", "#phi vs #eta",
140 20, -1, 1, 25, 0, 360,
141 "#eta", "#phi", "#Delta N/N",
142 AliMonitorHisto::kNormNone);
147 //_____________________________________________________________________________
148 void AliMonitorTPC::CreateBranches(TTree* tree)
150 // create a branch with TPC variables
152 tree->Branch("TPC", "AliMonitorDataTPC", &fData);
156 //_____________________________________________________________________________
157 void AliMonitorTPC::FillHistos(AliRunLoader* runLoader,
158 AliRawReader* rawReader, AliESD* esd)
160 // fill the TPC monitor histogrms
163 AliTPCRawStream input(rawReader);
164 while (input.Next()) {
165 fPadsCharge->Fill(input.GetSignal());
168 AliLoader* tpcLoader = runLoader->GetLoader("TPCLoader");
169 if (!tpcLoader) return;
172 tpcLoader->LoadRecPoints();
173 TTree* clusters = tpcLoader->TreeR();
174 if (!clusters) return;
175 AliTPCClustersRow* clustersRow = new AliTPCClustersRow;
176 clustersRow->SetClass("AliTPCclusterMI");
177 clusters->SetBranchAddress("Segment", &clustersRow);
179 for (Int_t i = 0; i < clusters->GetEntries(); i++) {
180 clusters->GetEntry(i);
182 fParam->AdjustSectorRow(clustersRow->GetID(), iSector, iRow);
183 if (iSector >= fParam->GetNInnerSector()) {
184 iRow += fParam->GetNRowLow();
185 iSector -= fParam->GetNInnerSector();
188 TObjArray* array = clustersRow->GetArray();
189 for (Int_t j = 0; j < array->GetEntriesFast(); j++) {
190 AliTPCclusterMI* cluster = (AliTPCclusterMI*) array->At(j);
191 fClustersCharge->Fill(cluster->GetQ());
192 fNClustersVsRow->Fill(iRow);
193 fNClustersVsSector->Fill(iSector);
196 fNClustersVsSector->ScaleErrorBy(10.);
199 tpcLoader->UnloadRecPoints();
203 for (Int_t i = 0; i < esd->GetNumberOfTracks(); i++) {
204 AliESDtrack* track = esd->GetTrack(i);
205 if (!track || ((track->GetStatus() | AliESDtrack::kTPCin) == 0)) continue;
209 track->GetInnerPxPyPz(pxyz);
210 TVector3 pTrack(pxyz);
211 Double_t p = pTrack.Mag();
212 Double_t pt = pTrack.Pt();
213 Double_t eta = pTrack.Eta();
214 Double_t phi = pTrack.Phi() * TMath::RadToDeg();
217 fTrackEta->Fill(eta);
218 fTrackPhi->Fill(phi);
220 fTrackEtaVsPhi->Fill(eta, phi);
221 fPtEtaVsPhi->Fill(eta, phi, pTrack.Pt());
223 fTrackNCl->Fill(track->GetTPCclusters(NULL));
224 fTrackDEdxVsP->Fill(p, track->GetTPCsignal());
226 fTrackDEdx->Fill(track->GetTPCsignal());
228 fData->SetData(i, pt, eta, phi);
230 fNTracks->Fill(nTracks);
231 fData->SetNTracks(nTracks);