Some warning going to error and viceversa:
[u/mrichter/AliRoot.git] / STEER / STEER / AliTransportMonitor.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 // Class that can be plugged in the simulation to monitor transport timing per 
17 // particle for each geometry volume.
18 //
19 //  andrei.gheata@cern.ch 
20
21 #include "AliTransportMonitor.h"
22
23 #include "TDatabasePDG.h"
24 #include "TCanvas.h"
25 #include "TMath.h"
26 #include "TFile.h"
27 #include "TROOT.h"
28 #include "THashList.h"
29 #include "TH2F.h"
30 #include "TGeoManager.h"
31 #include "AliPDG.h"
32
33 ClassImp(AliTransportMonitor)
34 ClassImp(AliTransportMonitor::AliTransportMonitorVol)
35 ClassImp(AliTransportMonitor::AliTransportMonitorVol::AliPMonData)
36
37 typedef AliTransportMonitor::AliTransportMonitorVol::AliPMonData PMonData;
38
39 //______________________________________________________________________________
40 AliTransportMonitor::AliTransportMonitorVol::AliTransportMonitorVol()
41                     :TNamed(),
42                      fNtypes(0),
43                      fTotalTime(0),
44                      fPData(0),
45                      fTimeRZ(0),
46                      fParticles()
47 {
48 // Default constructor
49 }
50
51 //______________________________________________________________________________
52 AliTransportMonitor::AliTransportMonitorVol::~AliTransportMonitorVol()
53 {
54 // Destructor
55   delete [] fPData;
56   delete fTimeRZ;
57 }
58
59 //______________________________________________________________________________
60 void AliTransportMonitor::AliTransportMonitorVol::StepInfo(
61                                     Int_t pdg,
62                                     Double_t energy, 
63                                     Double_t dt,
64                                     Double_t x, Double_t y, Double_t z)
65 {
66 // This method is called at each N steps to store timing info.
67   PMonData &data = GetPMonData(pdg);
68   data.fEdt += energy*dt;
69   data.fTime += dt;
70   fTotalTime += dt;
71   if (!fTimeRZ) {
72     Bool_t status = TH1::AddDirectoryStatus();
73     TH1::AddDirectory(kFALSE);
74     fTimeRZ = new TH2F("h2rz", "", 100, -5000., 5000, 100, 0., 1000.);
75     TH1::AddDirectory(status);
76   }
77   Double_t r = TMath::Sqrt(x*x+y*y);
78   fTimeRZ->Fill(z,r,dt);
79 }
80
81 //______________________________________________________________________________
82 PMonData &AliTransportMonitor::AliTransportMonitorVol::GetPMonData(Int_t pdg)
83 {
84 // Retrieve stored monitoring object for a given pdg type. If not existing 
85 // create one.
86
87   // The object could have been retrieved from file, in which case we have to 
88   // build the map.
89   PMonData *data;
90   if (fNtypes) {
91     if (fParticles.empty()) {
92       for (Int_t i=0; i<fNtypes; i++) {
93         data = &fPData[i];
94         fParticles[i] = data->fPDG;
95       }
96     }
97     ParticleMapIt_t it = fParticles.find(pdg);
98     if (it == fParticles.end()) {
99       data = &fPData[fNtypes]; 
100       data->fPDG = pdg;
101       fParticles[pdg] = fNtypes++;
102     } else {
103       data = &fPData[it->second];
104     }  
105   } else {
106      TDatabasePDG *pdgDB = TDatabasePDG::Instance();
107      if (!pdgDB->ParticleList()) AliPDG::AddParticlesToPdgDataBase();
108      Int_t size = pdgDB->ParticleList()->GetSize();
109      fPData = new PMonData[size];
110      data = &fPData[fNtypes];
111      data->fPDG = pdg;
112      fParticles[pdg] = fNtypes++;
113   }
114   return *data;
115 }
116
117 ClassImp(AliTransportMonitor)
118
119 //______________________________________________________________________________
120 AliTransportMonitor::AliTransportMonitor()
121                     :TObject(),
122                      fTotalTime(0),
123                      fTimer(),
124                      fVolumeMon(0)
125 {
126 // Default constructor
127 }
128
129 //______________________________________________________________________________
130 AliTransportMonitor::AliTransportMonitor(Int_t nvolumes)
131                     :TObject(),
132                      fTotalTime(0),
133                      fTimer(),
134                      fVolumeMon(0)
135 {
136 // Default constructor
137   fVolumeMon = new TObjArray(nvolumes);
138   fVolumeMon->SetOwner();
139   for (Int_t i=0; i<nvolumes; i++) {
140     AliTransportMonitorVol *volMon = new AliTransportMonitorVol();
141     if (gGeoManager) volMon->SetName(gGeoManager->GetListOfUVolumes()->At(i)->GetName());
142     fVolumeMon->Add(volMon);
143   }   
144 }
145
146 //______________________________________________________________________________
147 AliTransportMonitor::~AliTransportMonitor()
148 {
149 // Destructor
150   delete fVolumeMon;
151 }
152
153 //______________________________________________________________________________
154 void AliTransportMonitor::Print(Option_t *volName) const
155 {
156 // Inspect the timing statistics for a single volume or for all the setup
157   Int_t uid = -1;
158   Int_t ntotal = 0;
159   if (!fVolumeMon || !(ntotal=fVolumeMon->GetEntriesFast())) {
160      Info("Inspect", "Transport monitor is empty !");
161      return;
162   }   
163   if (strlen(volName)) {
164     TString svname = volName;
165     Int_t i = 0;
166     for (i=1; i<ntotal; i++) if (svname == fVolumeMon->At(i)->GetName()) break;
167     if (i==ntotal) {
168        Error("Inspect", "No monitoring info stored for volume %s", volName);
169        return;
170     }
171     uid = i;
172     AliTransportMonitorVol *volMon = (AliTransportMonitorVol*)fVolumeMon->At(uid);
173     Int_t ntypes = volMon->GetNtypes();
174     if (!ntypes) {
175       Info("Inspect", "No particles crossed volume %s", volName);
176       return;
177     }  
178     Double_t *timeperpart = new Double_t[ntypes];
179     Int_t *isort = new Int_t[ntypes];
180     Double_t timepervol = 0.;
181     for (i=0; i<ntypes; i++) {
182       timeperpart[i] = volMon->GetTime(i); 
183       timepervol += timeperpart[i];
184     }
185     printf("Volume %s: Transport time: %g%% of %g [s]\n", volMon->GetName(), 100.*timepervol/fTotalTime, fTotalTime);
186     TMath::Sort(ntypes, timeperpart, isort, kTRUE);
187     TString particle;
188     TDatabasePDG *pdgDB = TDatabasePDG::Instance();    
189     if (!pdgDB->ParticleList()) AliPDG::AddParticlesToPdgDataBase();
190     for (i=0; i<ntypes; i++)  {
191        timeperpart[i] /=  timepervol;
192        particle = pdgDB->GetParticle(volMon->GetPDG(isort[i]))->GetName();
193        printf("   %s: %g%%  mean energy: %g\n", particle.Data(), 100.*timeperpart[i], volMon->GetEmed(i));
194     }
195     if (volMon->GetHistogram()) {
196       TCanvas *c1 = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("crz");
197       if (!c1) c1 = new TCanvas("crz");
198       c1->cd();
199       volMon->GetHistogram()->GetXaxis()->SetTitle("z [cm]");
200       volMon->GetHistogram()->GetYaxis()->SetTitle("r [cm]");
201       volMon->GetHistogram()->SetTitle(Form("RZ plot weighted by time spent in %s",volMon->GetName()));
202       volMon->GetHistogram()->Draw();
203     }  
204     return;
205   }
206   // General view
207   TIter next(fVolumeMon);
208   AliTransportMonitorVol *volMon;
209   Int_t ncrossed = 0;
210   TH1F *hnames = new TH1F("volume timing", "relative volume timing", 3,0,3);
211   hnames->SetStats(0);
212   hnames->SetFillColor(38);
213   hnames->SetBit(TH1::kCanRebin);
214   while ((volMon=(AliTransportMonitorVol*)next())) {
215     if (volMon->GetNtypes()) {
216       hnames->Fill(volMon->GetName(), volMon->GetTotalTime());
217       ncrossed++;
218     }
219   }
220   
221   hnames->LabelsDeflate();
222   hnames->GetXaxis()->LabelsOption(">");
223   
224   TCanvas *c = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("cvol_timing");
225   if (!c) c = new TCanvas("cvol_timing");
226   c->cd();
227   c->SetLogy();
228   hnames->Draw();
229   
230   printf("=============================================================================\n");
231   printf("Effective transport time:  %6.2f minutes\n", fTotalTime/60.);
232   printf("Number of crossed volumes: %d from %d\n", ncrossed, fVolumeMon->GetEntriesFast());
233   printf("=============================================================================\n");  
234 }     
235
236 //______________________________________________________________________________
237 void AliTransportMonitor::DummyStep()
238 {
239 // Reset timer for zero-length steps
240    fTimer.Stop();
241    fTimer.Start(kTRUE);
242 }   
243
244 //______________________________________________________________________________
245 void AliTransportMonitor::StepInfo( Int_t volId,
246                                     Int_t pdg,
247                                     Double_t energy, 
248                                     Double_t x, Double_t y, Double_t z)
249 {
250 // This method is called at each N steps to store timing info.
251   fTimer.Stop();
252   Double_t dt = fTimer.RealTime();
253   fTotalTime += dt;
254   AliTransportMonitorVol *volMon = (AliTransportMonitorVol*)fVolumeMon->At(volId);
255   volMon->StepInfo(pdg,energy,dt,x,y,z);
256   fTimer.Start(kTRUE);
257 }
258
259 //______________________________________________________________________________
260 void AliTransportMonitor::Start()
261 {
262 // Start collecting timing information
263   if (fTotalTime > 0) {
264     Info("Start", "Cannot start twice");
265     return;
266   }
267   fTimer.Start(kTRUE);
268 }
269
270 //______________________________________________________________________________
271 void AliTransportMonitor::Stop()
272 {
273 // Stop the global timer
274   fTimer.Stop();
275 }
276
277 //______________________________________________________________________________
278 void AliTransportMonitor::Export(const char *fname)
279 {
280 // Export information to file
281   TFile *file = TFile::Open(fname, "RECREATE");
282   Write();
283   file->Write();
284   file->Close();
285 }  
286
287 //______________________________________________________________________________
288 AliTransportMonitor *AliTransportMonitor::Import(const char *fname)
289 {
290 // Import information from a file
291   TFile *file = TFile::Open(fname);
292   if (!file) {
293     ::Error("Import", "File %s could not be opened", fname);
294     return 0;
295   }
296   AliTransportMonitor *mon = (AliTransportMonitor *)file->Get("AliTransportMonitor");
297   if (!mon) {
298     ::Error("Import", "No AliTransportMonitor object found n file %s", fname);
299     return 0;
300   }
301   AliPDG::AddParticlesToPdgDataBase();
302   return mon;
303 }
304
305