Added static AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface function that operates...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALWsuCosmicRaySetUp.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2010, 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 //  Wsu Cosmic Ray SetUp                                                     //
21 //  This class contains the description of the  Wsu Cosmic Ray SetUp         //
22 //  external volume                                                          //
23 //                                                                           //
24 //Begin_Html
25 /*
26 <img src="picts/AliEMCALWsuCosmicRaySetUpClass.gif">
27 </pre>
28 <br clear=left>
29 <font size=+2 color=red>
30 <p>The responsible person for this module is
31 <a href="mailto:pavlinov@physics.wayne.edu">Alexei Pavlino, WSU</a>.
32 </font>
33 <pre>
34 */
35 //End_Html
36 //                                                                           //
37 //                                                                           //
38 ///////////////////////////////////////////////////////////////////////////////
39
40 #include <TVirtualMC.h>
41 #include <TROOT.h>
42 #include <TDatabasePDG.h>
43 #include <TBrowser.h>
44 #include <TLorentzVector.h>
45 #include <TParticle.h>
46 #include <TList.h>
47 #include <TH1.h>
48 #include <TH2.h>
49 #include <TAxis.h>
50
51 #include "AliEMCALWsuCosmicRaySetUp.h"
52 //#include "AliMagF.h"
53 #include "AliStack.h"
54 #include "AliRun.h"
55 #include "AliMC.h"
56
57 #include "AliEMCALHistoUtilities.h"
58
59 using namespace std;
60
61 typedef AliEMCALHistoUtilities hist;
62
63 TDatabasePDG *pdg = 0; 
64
65 Int_t gid=0;
66
67 ClassImp(AliEMCALWsuCosmicRaySetUp)
68
69 //  TList *fLHists=0, *ll=0; 
70
71 //_____________________________________________________________________________
72   AliEMCALWsuCosmicRaySetUp::AliEMCALWsuCosmicRaySetUp(): AliModule()
73                                                       //,fMasterVolume()
74                                                     ,fLHists(0),fMasterVolume()
75 {
76   //
77   // Default constructor
78   //
79 }
80  
81 //_____________________________________________________________________________
82 AliEMCALWsuCosmicRaySetUp::AliEMCALWsuCosmicRaySetUp(const char *name, const char *title)
83        : AliModule(name,title)
84                                                     //,fMasterVolume()
85                                                     ,fLHists(0),fMasterVolume()
86 {
87   //
88   // Standard constructor of the  Wsu Cosmic Ray SetUp external volume
89   //
90   //PH  SetMarkerColor(7);
91   //PH  SetMarkerStyle(2);
92   //PH  SetMarkerSize(0.4);
93 }
94  
95 //_____________________________________________________________________________
96 void AliEMCALWsuCosmicRaySetUp::CreateGeometry()
97 {
98   //
99   // Create the geometry of the Alice external body
100   //
101   //Begin_Html
102   /*
103     <img src="picts/AliEMCALWsuCosmicRaySetUpTree.gif">
104   */
105   //End_Html
106   pdg = TDatabasePDG::Instance();
107
108   // Master Volume
109   fMasterVolume[0] = fMasterVolume[1] = 35.0;
110   fMasterVolume[2] = 450.;
111
112   Int_t *idtmed = fIdtmed->GetArray()+1;
113   int idAir = idtmed[0];
114   gMC->Gsvolu(GetName(),"BOX",idAir, fMasterVolume,3); // Master volume
115   //
116   // Sc counters
117   //
118   Float_t sc[3]; // tube
119   sc[0] = 0.0;
120   sc[1] = 10.0;
121   sc[2] = 1.0; // thicness of Sc is 2 cm
122   Float_t zsc[3] = {10.,330.6, 810.1}; 
123   int idSC = idtmed[1];
124   gMC->Gsvolu("SCOU","TUBE",idSC, sc,3); // Master volume
125   printf(" idtmed[0] %i idtmed[1] %i \n", idtmed[0] , idtmed[1]); 
126   Int_t idRot=0; // no rotation
127   for(Int_t i=0; i<3; i++) {
128     Float_t zpos = zsc[i] - fMasterVolume[2];
129     gMC->Gspos("SCOU", i+1, "WSUC", 0.0, 0.0, zpos, idRot, "ONLY"); 
130   }
131   //
132   // Dead end : Dec 2,2010
133   //
134   Float_t zbox[3]={30., 30.0, 0.1};
135   gMC->Gsvolu("SEND","BOX",idAir, zbox,3); // Master volume
136   gMC->Gspos("SEND", 1, "WSUC", 0.0, 0.0, 448.0, idRot, "ONLY"); 
137   // Hists
138   fLHists = new TList;
139   fLHists->SetName("hists");
140   //
141   //AliMC *ALIMC  = dynamic_cast<AliMC *>(gMC);
142   //AliGenBox* gB = dynamic_cast<AliGenBox *>(ALIMC->Generator());
143   //Double_t p = gB->
144   Double_t pmom=1.5; 
145   fLHists->Add(BookKineHists(pmom,"primeKineHists"));
146   fLHists->Add(BookKineHists(pmom,"endKineHists"));
147   fLHists->Add(BookKineHists(pmom,"secondaryKineHists"));
148   //ll = BookKineHists(1.,"kineHists");
149   //gROOT->GetListOfBrowsables()->Add(ll);
150 }
151  
152 //_____________________________________________________________________________
153 void AliEMCALWsuCosmicRaySetUp::CreateMaterials()
154 {
155 // Create materials and media
156   Int_t   isxfld = 0;
157   Float_t sxmgmx = 0., deemax = 0.1;  
158   // AIR
159   Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
160   Float_t zAir[4]={6.,7.,8.,18.};
161   Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
162   Float_t dAir = 1.20479E-3;
163   AliMixture(1,"Air     $",aAir,zAir,dAir,4,wAir);
164
165   // --- The polysterene scintillator (CH) ---
166   Float_t aP[2] = {12.011, 1.00794} ;
167   Float_t zP[2] = {6.0, 1.0} ;
168   Float_t wP[2] = {1.0, 1.0} ;
169   Float_t dP = 1.032 ;
170   AliMixture(2, "Polystyrene$", aP, zP, dP, -2, wP) ;
171
172   //
173   AliMedium(1,"Air     $",1,0,isxfld,sxmgmx,10,-1,-0.1,0.1 ,-10);
174   AliMedium(2, "Scintillator$", 2, 1,
175             isxfld, sxmgmx, 10.0, 0.001, deemax, 0.001, 0.001, 0, 0) ;
176   // cuts
177   DefineCuts(1);
178   DefineCuts(2);
179 }
180  
181 void AliEMCALWsuCosmicRaySetUp::DefineCuts(const Int_t idtmed)
182 {
183   // Dec 2,2010 : it works
184   Float_t cutgam=10.e-5; // 100 kev;
185   Float_t cutele=10.e-5; // 100 kev;
186   gMC->Gstpar(idtmed,"CUTGAM", cutgam);
187   gMC->Gstpar(idtmed,"CUTELE", cutele); // 1MEV -> 0.1MEV; 15-aug-05
188   gMC->Gstpar(idtmed,"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
189   gMC->Gstpar(idtmed,"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
190   // --- Generate explicitly delta rays in Lead ---
191   gMC->Gstpar(idtmed, "LOSS", 3) ;
192   gMC->Gstpar(idtmed, "DRAY", 1) ;
193   gMC->Gstpar(idtmed, "DCUTE", cutele) ;
194   gMC->Gstpar(idtmed, "DCUTM", cutele) ;
195 }
196  
197 void AliEMCALWsuCosmicRaySetUp::StepManager(void)
198 {
199   // Dec 1,2010
200   static Int_t pri=1;
201   static TString curVolName="";
202   static TLorentzVector pos;  // Lorentz vector of the track current position.
203   static TLorentzVector mom;  // Lorentz vector of the track current momentum.
204
205   if(pri>=2) printf("<I> AliEMCALWsuCosmicRaySetUp::StepManager %s \n", gMC->CurrentVolName());
206   Int_t tracknumber =  gAlice->GetMCApp()->GetCurrentTrackNumber();
207   //  Int_t parent=0;
208   TParticle* part=0;
209   curVolName = gMC->CurrentVolName();
210   if(curVolName.Contains("SEND")) {
211     gMC->TrackMomentum(mom);
212     gMC->TrackPosition(pos);
213     if(pri>=2) printf(" %s tracknumber %i p %f \n", curVolName.Data(), tracknumber, mom.P());
214     if(pri>=2) printf(" x %f y %f z %f \n", pos[0], pos[1], pos[2]);
215     if(gMC->IsTrackEntering()) { // primary only TList *l = GetLhists(1);
216       TList *l = 0;
217       if(tracknumber==0){
218         l = GetLhists(1);
219         part=gAlice->GetMCApp()->Particle(0);
220         gid = pdg->ConvertPdgToGeant3(part->GetPdgCode()); 
221         hist::FillH1(l, 1, double(gid));
222       } else {
223         l = GetLhists(2);
224         part=gAlice->GetMCApp()->Particle(tracknumber);
225         gid = pdg->ConvertPdgToGeant3(part->GetPdgCode()); 
226         hist::FillH1(l, 1, double(gid));
227       }
228       Int_t ic = 2;
229       hist::FillH1(l, ic++, mom.P());
230       hist::FillH1(l, ic++, mom.Eta());
231       hist::FillH1(l, ic++, TVector2::Phi_0_2pi(mom.Phi())*TMath::RadToDeg() );
232       hist::FillH1(l, ic++, mom.Theta()*TMath::RadToDeg());
233     }
234   }
235 }
236
237 void AliEMCALWsuCosmicRaySetUp::FinishEvent()
238 {
239   // Dec 2,2010
240   int ic=0;
241
242   TList *l = GetLhists(0);
243   //  TList *l = ll;
244   AliStack* st =  AliRunLoader::Instance()->Stack();
245   TParticle *p = st->Particle(0);
246   gid = pdg->ConvertPdgToGeant3(p->GetPdgCode());
247
248   ic=1;
249   hist::FillH1(l, ic++, double(gid));
250   hist::FillH1(l, ic++, p->P());
251   hist::FillH1(l, ic++, p->Eta());
252   hist::FillH1(l, ic++, TVector2::Phi_0_2pi(p->Phi())*TMath::RadToDeg() );
253   hist::FillH1(l, ic++, p->Theta()*TMath::RadToDeg());
254 }
255
256 TList* AliEMCALWsuCosmicRaySetUp::BookKineHists(const Double_t p , const Char_t *opt)
257 {
258   // Dec 2,2010
259   gROOT->cd();
260   TH1::AddDirectory(1);
261
262   TH1 *h = 0, *hgid=0;
263   Int_t nphi=180, nmax=1100;
264   Double_t phimin=0.0, phimax=360.;
265   Double_t pmax=110.;
266   if(p>0.1) pmax = 1.1*p;
267   h = new TH1F("00_hNPrim"," number of primary particles ", 10, 0.5, 10.5);
268   hgid = new TH1F("01_hGidprim","Geant Id of primary particles", 16, 0.5, 16.5);
269   new TH1F("02_hPmomPrim","p of primary particles", nmax, 0.0, pmax);
270   new TH1F("03_hEtaPrim","#eta primary particles", 80, 0.0, 8.0);
271   new TH1F("04_hPhiPrim","#phi primary particles", nphi,phimin,phimax);
272   new TH1F("05_hThetaPrim","#theta primary particles", 90, 0.0, 90.);
273   //
274   TAxis *xax=hgid->GetXaxis();
275   xax->SetBinLabel(1,"#gamma");
276   xax->SetBinLabel(2,"e^{+}");
277   xax->SetBinLabel(3,"e^{-}");
278   xax->SetBinLabel(4,"#nu");
279   xax->SetBinLabel(5,"#mu^{+}");
280   xax->SetBinLabel(6,"#mu^{-}");
281   xax->SetBinLabel(7,"#pi^{0}");
282   xax->SetBinLabel(8,"#pi^{+}");
283   xax->SetBinLabel(9,"#pi^{-}");
284   xax->SetBinLabel(10,"K_{L}");
285   xax->SetBinLabel(11,"K^{+}");
286   xax->SetBinLabel(12,"K^{-}");
287   xax->SetBinLabel(13,"n");
288   xax->SetBinLabel(14,"p");
289   xax->SetBinLabel(15,"#bar{p}");
290   xax->SetBinLabel(16,"K_{s}");
291   //  hgid->SetBinLabel(,"");
292
293   TList *l = 0;
294   l = hist::MoveHistsToList(opt, 1);
295   if(strlen(opt)) hist::AddToNameAndTitleToList(l, opt, opt);
296
297   TH1::AddDirectory(0);
298
299   return l;
300 }
301
302 void AliEMCALWsuCosmicRaySetUp::Browse(TBrowser* b)
303 {
304   if(fLHists) b->Add(fLHists);
305   TObject::Browse(b);
306 }
307 /*
308 */