ITS/UPGRADE/readHit.C

   1 void readHit(int nev=-1,int evStart=0)
   2 {
   3   gROOT->SetStyle("Plain");
   4
   5   gSystem->Load("libITSUpgradeBase");
   6   gSystem->Load("libITSUpgradeSim");
   7
   8   TH2F *xyGlob = new TH2F("xyGlob"," X - Y Global coordinates ",100,-50,50,100,-50,50);
   9   xyGlob->SetXTitle("cm");
  10   xyGlob->SetMarkerStyle(7);
  11   TH1F *zGlob  = new TH1F("zGlob", " Z Global coordinates ",200, -100,100 );
  12   zGlob->SetXTitle("cm");
  13
  14   Int_t nbins=100;
  15   Double_t xmin=0;
  16   Double_t xmax=1e-04;
  17
  18   gAlice=NULL;
  19   AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
  20   runLoader->LoadgAlice();
  21
  22   gAlice = runLoader->GetAliRun();
  23
  24   runLoader->LoadHeader();
  25   runLoader->LoadKinematics();
  26   runLoader->LoadSDigits();
  27   runLoader->LoadHits();
  28
  29   AliGeomManager::LoadGeometry("geometry.root");
  30   AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
  31
  32   AliLoader *dl = runLoader->GetDetectorLoader("ITS");
  33
  34
  35
  36   Int_t nLayers = gm->GetNLayers();
  37
  38   TH1D **hDeLoss = new TH1D*[nLayers];
  39   for(Int_t i=0; i< nLayers; i++ ) {
  40     hDeLoss[i] = new TH1D(Form("hDeLossl%i",i),Form("E loss distribution [ Layer %i] ",i),nbins,xmin,xmax);
  41     hDeLoss[i]->SetXTitle("GeV");
  42   }
  43
  44   int nevTot = (Int_t)runLoader->GetNumberOfEvents();
  45   printf("N Events : %i \n",nevTot);
  46   evStart = evStart<nevTot ? evStart : nevTot-1;
  47   if (evStart<0) evStart = 0;
  48   //
  49   int lastEv = nev<0 ? nevTot : evStart+nev;
  50   if (lastEv > nevTot) lastEv = nevTot;
  51   //
  52
  53   //HIT INIT
  54
  55   TTree *hitTree = 0x0;
  56   TClonesArray *hitList=new TClonesArray("AliITSUHit");
  57
  58   for (Int_t iEvent = evStart; iEvent < lastEv; iEvent++) {
  59
  60     printf("\nEvent\t%d\n",iEvent);
  61
  62     runLoader->GetEvent(iEvent);
  63
  64     hitTree=dl->TreeH();
  65     hitTree->SetBranchAddress("ITS",&hitList);
  66     for(Int_t iEnt=0;iEnt<hitTree->GetEntries();iEnt++){//entries loop degli hits
  67       hitTree->GetEntry(iEnt);
  68       for(Int_t iHit=0; iHit<hitList->GetEntries();iHit++){
  69
  70         AliITSUHit *pHit = (AliITSUHit*)hitList->At(iHit);
  71         int id = pHit->GetModule();
  72         int lr = gm->GetLayer(id);
  73         int ld = gm->GetLadder(id);
  74         //
  75         //      if(pHit->GetParticle()->IsPrimary()){
  76         Double_t xg,yg,zg=0.;
  77         Double_t xg0,yg0,zg0=0.,tg0;
  78         pHit->GetPositionG(xg,yg,zg);
  79         pHit->GetPositionG0(xg0,yg0,zg0,tg0);
  80         xyGlob->Fill(xg,yg);
  81         zGlob->Fill(zg);
  82         printf("Module %5d | Lr:%2d Ladder: %3d, X:[%+.5e:%+.5e] Y:[%+.5e:%+.5e] Z:[%+.5e %+.5e] DE: %.5e TrackID: %d\n",id,lr,ld,
  83                xg0,xg,yg0,yg,zg0,zg,pHit->GetIonization(),pHit->GetTrack());
  84         hDeLoss[lr]->Fill(pHit->GetIonization());
  85         //      } // is primary
  86       }//loop hit
  87     }//entryloopHitList
  88
  89   }//event loop
  90
  91   TCanvas *xyCanv =  new TCanvas("xyCanv","Hit X-Y positions",500,500);
  92   xyCanv->cd();
  93   xyGlob->Draw();
  94   TCanvas *zCanv =  new TCanvas("zCanv","Hit Z positions",500,500);
  95   zCanv->cd();
  96   zGlob->Draw();
  97
  98   TCanvas *c = new TCanvas("c","E loss  distribution per layer",1000,800);
  99   c->Divide(nLayers&0x1 ?  (nLayers/2+1) : nLayers/2,2);
 100   for(Int_t ip =1; ip<=nLayers; ip++){
 101     c->cd(ip);
 102     hDeLoss[ip-1]->Draw();
 103   }
 104
 105 }