Updated SNM Glauber fit
[u/mrichter/AliRoot.git] / macros / printProfile.C
1 void printProfile(char* g3 = "geant3/timing_200_woZDC.root", char* g4 ="geant4/timing_200_woZDC.root")
2 {
3   // Create monitors 
4   AliTransportMonitor* mon3 = AliTransportMonitor::Import(g3);
5   AliTransportMonitor* mon4 = AliTransportMonitor::Import(g4);
6   //
7   // Create timing histos and volume lists
8   mon3->Print();
9   tH3 = (TH1F*) volume_timing->Clone();
10   TObjArray* volumes3 = mon3->GetVolumes(); 
11   mon4->Print();
12   tH4 = (TH1F*) volume_timing->Clone();
13   TObjArray* volumes4 = mon4->GetVolumes(); 
14   //
15   // Normalise the geant 3 histo
16   // Fraction of time in %
17   Float_t ri = tH3->Integral();
18   tH3->Scale(100./ri);
19   //
20   Float_t cumulant = 0.;  // cumulant
21   Int_t   i = 1;          // position in list
22   //
23   // print up to 90% of total time
24   //
25   printf(" Pos.               Volume  t (%) ct (%)  Time_G3  Time_G4 perStepG3 perStepG4   StepsG3   StepsG4    te+e-     tgam     thad     te+e-    tgam     thad\n"); 
26   while (cumulant < 90.) {
27     Float_t rt = tH3->GetBinContent(i);
28     cumulant += rt;
29     // extract the corresponding geant3 and geant4 volumes
30     char* volN = (tH3->GetXaxis())->GetBinLabel(i);
31   AliTransportMonitor::AliTransportMonitorVol* vol4 = (AliTransportMonitor::AliTransportMonitorVol*) 
32     volumes4->FindObject(volN);
33   AliTransportMonitor::AliTransportMonitorVol* vol3 = 0;
34     if (i > 1) {
35       vol3 = (AliTransportMonitor::AliTransportMonitorVol*) 
36         volumes3->FindObject(volN);
37     } else {
38       vol3 = (AliTransportMonitor::AliTransportMonitorVol*) 
39         volumes3->At(1);
40      }
41     TString s(volN);
42     TString s20(s(0,20));
43     // Total time
44     
45     Float_t t4[5]; 
46     t4[0] = vol4->GetTotalTime();
47     Float_t t3[5]; 
48     t3[0] = vol3->GetTotalTime();
49
50     // For e+/e-
51     t4[1] = ElectronTime(vol4);
52     t3[1] = ElectronTime(vol3);
53  
54     // For photons
55     t4[2] = PhotonTime(vol4);
56     t3[2] = PhotonTime(vol3);
57    
58     // Hadrons
59     t4[3] = HadronTime(vol4);
60     t3[3] = HadronTime(vol3);
61    
62     printf("%5d %20s %6.2f %6.2f %8.3f %8.3f %8.3e %8.3e %8.3e %8.3e %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n", 
63            i,
64            s20.Data(),
65            rt, cumulant, t3[0], t4[0], t3[0]/(1.+vol3->GetNSteps()), 
66            t4[0]/(1.+vol4->GetNSteps()), vol3->GetNSteps(), vol4->GetNSteps(),
67            t3[1]/t3[0] * 100., t3[2]/t3[0] * 100., t3[3]/t3[0] * 100.,
68            t4[1]/t4[0] * 100., t4[2]/t4[0] * 100., t4[3]/t4[0] * 100.
69
70            );
71     i++;
72   }
73 }
74
75 Float_t ElectronTime(AliTransportMonitor::AliTransportMonitorVol* vol)
76 {
77    Int_t nt = vol->GetNtypes();
78    Float_t t = 0.;
79    for (Int_t i = 0; i < nt; i++) {
80      Int_t pdg = TMath::Abs(vol->GetPDG(i));
81      if (pdg == 11) t += vol->GetTime(i);
82    }
83
84    return t;
85 }
86
87 Float_t PhotonTime(AliTransportMonitor::AliTransportMonitorVol* vol)
88 {
89    Int_t nt = vol->GetNtypes();
90    Float_t t = 0.;
91    for (Int_t i = 0; i < nt; i++) {
92      Int_t pdg = TMath::Abs(vol->GetPDG(i));
93      if (pdg == 22) t += vol->GetTime(i);
94    }
95    return t;
96 }
97
98 Float_t HadronTime(AliTransportMonitor::AliTransportMonitorVol* vol)
99 {
100    Int_t nt = vol->GetNtypes();
101    Float_t t = 0.;
102    for (Int_t i = 0; i < nt; i++) {
103      Int_t pdg = TMath::Abs(vol->GetPDG(i));
104      if (pdg == 211 || pdg == 2212 || pdg == 2112 || pdg == 130 || pdg == 321) t += vol->GetTime(i);
105    }
106    return t;
107 }
108
109