EMCAL/EMCALbase/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   TVirtualMC::GetMC()->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   TVirtualMC::GetMC()->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     TVirtualMC::GetMC()->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   TVirtualMC::GetMC()->Gsvolu("SEND","BOX",idAir, zbox,3); // Master volume
 136   TVirtualMC::GetMC()->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 *>(TVirtualMC::GetMC());
 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   TVirtualMC::GetMC()->Gstpar(idtmed,"CUTGAM", cutgam);
 187   TVirtualMC::GetMC()->Gstpar(idtmed,"CUTELE", cutele); // 1MEV -> 0.1MEV; 15-aug-05
 188   TVirtualMC::GetMC()->Gstpar(idtmed,"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
 189   TVirtualMC::GetMC()->Gstpar(idtmed,"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
 190   // --- Generate explicitly delta rays in Lead ---
 191   TVirtualMC::GetMC()->Gstpar(idtmed, "LOSS", 3) ;
 192   TVirtualMC::GetMC()->Gstpar(idtmed, "DRAY", 1) ;
 193   TVirtualMC::GetMC()->Gstpar(idtmed, "DCUTE", cutele) ;
 194   TVirtualMC::GetMC()->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", TVirtualMC::GetMC()->CurrentVolName());
 206   Int_t tracknumber =  gAlice->GetMCApp()->GetCurrentTrackNumber();
 207   //  Int_t parent=0;
 208   TParticle* part=0;
 209   curVolName = TVirtualMC::GetMC()->CurrentVolName();
 210   if(curVolName.Contains("SEND")) {
 211     TVirtualMC::GetMC()->TrackMomentum(mom);
 212     TVirtualMC::GetMC()->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(TVirtualMC::GetMC()->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 * 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   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 */