ZDC/ZDCtest.C

   1 void ZDCtest (Int_t detector=0, Int_t evTot = 0)
   2 {
   3    delete gAlice;
   4    gAlice=0;
   5 // Dynamically link some shared libs
   6    if (gClassTable->GetID("AliRun") < 0) {
   7       gROOT->LoadMacro("loadlibs.C");
   8       loadlibs();
   9    }
  10
  11 // Connect the Root Galice file containing Geometry, Kine, Hits and Digits
  12    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
  13     if (!file) {
  14         printf("\n Creating galice.root \n");
  15         file = new TFile("galice.root");
  16     } else {
  17         printf("\n galice.root found in file list");
  18     }
  19     file->ls();
  20
  21 // Get AliRun object from file or create it if not on file
  22    if (!gAlice) {
  23       gAlice = (AliRun*)file->Get("gAlice");
  24       if (gAlice) printf("AliRun object found on file\n");
  25       if (!gAlice) {
  26             printf("\n create new gAlice object");
  27             gAlice = new AliRun("gAlice","Alice test program");
  28         }
  29    }
  30    file->ls();
  31
  32 // Create some histograms
  33
  34    TH2F *hspotzn  = new TH2F("hspotzn","Y vs X on ZN front face",100,-4.,4.,100,-4.,4.);
  35    hspotzn -> SetXTitle("X on ZN");
  36    hspotzn -> SetYTitle("Y on ZN");
  37    TH2F *hspotzp  = new TH2F("hspotzp","Y vs X on ZP front face",100,-12.,12.,100,-12.,12.);
  38    hspotzp -> SetXTitle("X on ZP");
  39    hspotzp -> SetYTitle("Y on ZP");
  40    TH2F *hspotzem = new TH2F("hspotzem","Y vs X on ZEM front face",100,-4.,4.,100,-4.,4.);
  41    hspotzem -> SetXTitle("X on ZEM");
  42    hspotzem -> SetYTitle("Y on ZEM");
  43
  44    TH1F *hEzn    = new TH1F("hEzn","Energy deposited in ZN",100,0,4.e3);
  45    hEzn -> SetXTitle("E (GeV)");
  46    TH1F *hLzn    = new TH1F("hLzn", "Total light in ZN    ",100,0,4.e3);
  47    hLzn -> SetXTitle("phe");
  48    TH1F *hPMCzn  = new TH1F("hPMCzn", "Light in common PM    ",100,0,1.e3);
  49    hPMCzn -> SetXTitle("phe");
  50    TH1F *hPMQ1zn = new TH1F("hPMQ1zn","Light in quadrant 1 PM",100,0,1.e3);
  51    hPMQ1zn -> SetXTitle("phe");
  52    TH1F *hPMQ2zn = new TH1F("hPMQ2zn","Light in quadrant 2 PM",100,0,1.e3);
  53    hPMQ2zn -> SetXTitle("phe");
  54    TH1F *hPMQ3zn = new TH1F("hPMQ3zn","Light in quadrant 3 PM",100,0,1.e3);
  55    hPMQ3zn -> SetXTitle("phe");
  56    TH1F *hPMQ4zn = new TH1F("hPMQ4zn","Light in quadrant 4 PM",100,0,1.e3);
  57    hPMQ4zn -> SetXTitle("phe");
  58
  59    TH1F *hEzp    = new TH1F("hEzp","Energy deposited in ZP",100,0,4.e3);
  60    hEzp -> SetXTitle("E (GeV)");
  61    TH1F *hLzp    = new TH1F("hLzp", "Total light in ZP    ",100,0,4.e3);
  62    hLzp -> SetXTitle("phe");
  63    TH1F *hPMCzp  = new TH1F("hPMCzp","Light in common PM     ",100,0,1.e3);
  64    hPMCzp -> SetXTitle("phe");
  65    TH1F *hPMQ1zp = new TH1F("hPMQ1zp","Light in quadrant 1 PM",100,0,1.e3);
  66    hPMQ1zp -> SetXTitle("phe");
  67    TH1F *hPMQ2zp = new TH1F("hPMQ2zp","Light in quadrant 2 PM",100,0,1.e3);
  68    hPMQ2zp -> SetXTitle("phe");
  69    TH1F *hPMQ3zp = new TH1F("hPMQ3zp","Light in quadrant 3 PM",100,0,1.e3);
  70    hPMQ3zp -> SetXTitle("phe");
  71    TH1F *hPMQ4zp = new TH1F("hPMQ4zp","Light in quadrant 4 PM",100,0,1.e3);
  72    hPMQ4zp -> SetXTitle("phe");
  73
  74
  75    TH1F *hEzem   = new TH1F("hEzem","Energy deposited in ZEM",100,0,5.e3);
  76    hEzem -> SetXTitle("E (GeV)");
  77    TH1F *hPMzem  = new TH1F("hPMzem","Light produced in ZEM PM",100,0,5.e3);
  78    hPMzem -> SetXTitle("phe");
  79
  80    //
  81
  82    TH1F *dPMCzn  = new TH1F("dPMCzn","Common PM     ",100,0,1.e3);
  83    dPMCzn -> SetXTitle("ADC channels");
  84    TH1F *dPMQ1zn = new TH1F("dPMQ1zn","Quadrant 1 PM",100,0,1.e3);
  85    dPMQ1zn -> SetXTitle("ADC channels");
  86    TH1F *dPMQ2zn = new TH1F("dPMQ2zn","Quadrant 2 PM",100,0,1.e3);
  87    dPMQ2zn -> SetXTitle("ADC channels");
  88    TH1F *dPMQ3zn = new TH1F("dPMQ3zn","Quadrant 3 PM",100,0,1.e3);
  89    dPMQ3zn -> SetXTitle("ADC channels");
  90    TH1F *dPMQ4zn = new TH1F("dPMQ4zn","Quadrant 4 PM",100,0,1.e3);
  91    dPMQ4zn -> SetXTitle("ADC channels");
  92    TH1F *dZN     = new TH1F("dZN","Total light in ZN",100,0,1.e3);
  93    dZN -> SetXTitle("ADC channels");
  94
  95    TH1F *dPMCzp  = new TH1F("dPMCzp","Common PM     ",100,0,1.e3);
  96    dPMCzp -> SetXTitle("ADC channels");
  97    TH1F *dPMQ1zp = new TH1F("dPMQ1zp","Quadrant 1 PM",100,0,1.e3);
  98    dPMQ1zp -> SetXTitle("ADC channels");
  99    TH1F *dPMQ2zp = new TH1F("dPMQ2zp","Quadrant 2 PM",100,0,1.e3);
 100    dPMQ2zp -> SetXTitle("ADC channels");
 101    TH1F *dPMQ3zp = new TH1F("dPMQ3zp","Quadrant 3 PM",100,0,1.e3);
 102    dPMQ3zp -> SetXTitle("ADC channels");
 103    TH1F *dPMQ4zp = new TH1F("dPMQ4zp","Quadrant 4 PM",100,0,1.e3);
 104    dPMQ4zp -> SetXTitle("ADC channels");
 105    TH1F *dZP     = new TH1F("dZP","Total light in ZP",100,0,1.e3);
 106    dZP -> SetXTitle("ADC channels");
 107
 108    TH1F *dZEM     = new TH1F("dZEM","Total light in ZEM",100,0,1.e3);
 109    dZEM -> SetXTitle("ADC channels");
 110
 111 //
 112 //   Loop over events
 113 //
 114 // NB -> Il TClonesArray delle particelle va inizializzato prima del loop
 115    TParticle   *particle;
 116    TClonesArray *Particles = gAlice->Particles();
 117
 118    for (Int_t evNumber=0; evNumber<evTot; evNumber++){
 119
 120 // Import the Kine and Hits Trees for the event evNumber in the file
 121    Int_t nparticles = gAlice->GetEvent(evNumber);
 122    if (nparticles <= 0) return;
 123    printf("\n   --- nparticles = %d\n",nparticles);
 124
 125    Float_t energy, EtotZN=0, EtotZP=0, LightCzn=0, LightCzp=0, LtotZN=0, LtotZP=0;
 126    Int_t nbytes=0, nbytesd=0, ipart, nhits, ndigits, pdgcode, ADCzn, ADCzp;
 127    AliZDCHit   *ZDChit;
 128    AliZDCDigit *ZDCdigit;
 129
 130 // Get pointers to Alice detectors and Hits containers
 131    AliDetector *ZDC  = gAlice->GetDetector("ZDC");
 132    if (ZDC) {
 133      TClonesArray *ZDChits    = ZDC->Hits();
 134 //     TClonesArray *ZDCdigits  = ZDC->Digits();
 135    }
 136
 137 // # of entries in Hits tree
 138    TTree *TH = gAlice->TreeH();
 139    Int_t ntracks = TH->GetEntries();
 140
 141 // # of entries in Digits tree
 142 //   TTree *TD = gAlice->TreeD();
 143 //   Int_t ndigen = TD->GetEntries();
 144
 145 //   gAlice->ResetDigits();
 146 //   nbytesd += TD->GetEvent(ndigen-1);
 147
 148 //   if (ZDC) {
 149 //     ndigits = ZDCdigits->GetEntries();
 150 //     printf("\n       Digits Tree --- # of entries: %d; # of digits: %d\n",ndigen, ndigits);
 151 //   }
 152 //   for(Int_t digit=0; digit<ndigits; digit++) {
 153 //     ZDCdigit = (AliZDCDigit*)ZDCdigits->UncheckedAt(digit);
 154 //     printf("\n       Digit# %d, fDetector = %d, fVolume = %d, fADCValue = %f\n",
 155 //            digit,ZDCdigit->fDetector,ZDCdigit->fQuadrant,ZDCdigit->fADCValue);
 156 //   }
 157
 158 // Start loop on tracks in the hits containers
 159       for (Int_t track=0; track<ntracks; track++) {
 160          gAlice->ResetHits();
 161          nbytes += TH->GetEvent(track);
 162
 163          if (ZDC) {
 164 //          nhits = ZDChits->GetEntries();
 165 //          nhits = ZDChits->GetLast()+1;
 166            nhits   = ZDChits->GetEntriesFast();
 167 //           printf("\n Hits Tree --- Event %d track %d nhits %d\n",evNumber,track,nhits);
 168
 169            particle = (TParticle*)Particles->UncheckedAt(track);
 170 //         pdgcode = particle->GetPdgCode();
 171 //         printf("\nParticle %d\n",pdgcode);
 172
 173            for(Int_t hit=0; hit<nhits; hit++) {
 174              ZDChit   = (AliZDCHit*)ZDChits->UncheckedAt(hit);
 175
 176              // Print of the hits
 177 //           printf("\nHit # %d, fVolume = %d %d\n",hit,ZDChit->fVolume[0],ZDChit->fVolume[1]);
 178 //           printf("Primary energy = %f, Secondary Flag = %d\n",ZDChit->fPrimKinEn,ZDChit->fSFlag);
 179 //           printf("Impact point -> %f %f\n",ZDChit->fXImpact,ZDChit->fYImpact);
 180 //           printf("Energy = %f, Light in quadrant = %f, Light in common PM = %f\n\n",
 181 //                   ZDChit->fEnergy,ZDChit->fLightPMQ,ZDChit->fLightPMC);
 182
 183              // Filling histos
 184              if(ZDChit->fVolume[0]==1) { //ZN
 185                if(ZDChit->fVolume[1]==1)hPMQ1zn->Fill(ZDChit->fLightPMQ);
 186                if(ZDChit->fVolume[1]==2)hPMQ2zn->Fill(ZDChit->fLightPMQ);
 187                if(ZDChit->fVolume[1]==3)hPMQ3zn->Fill(ZDChit->fLightPMQ);
 188                if(ZDChit->fVolume[1]==4)hPMQ4zn->Fill(ZDChit->fLightPMQ);
 189                EtotZN   += ZDChit->fEnergy;
 190                LtotZN   += (ZDChit->fLightPMQ) + (ZDChit->fLightPMC);
 191                LightCzn += ZDChit->fLightPMC;
 192                hspotzn->Fill(ZDChit->fXImpact,ZDChit->fYImpact);
 193              }
 194              if(ZDChit->fVolume[0]==2) { //ZP
 195                if(ZDChit->fVolume[1]==1)hPMQ1zp->Fill(ZDChit->fLightPMQ);
 196                if(ZDChit->fVolume[1]==2)hPMQ2zp->Fill(ZDChit->fLightPMQ);
 197                if(ZDChit->fVolume[1]==3)hPMQ3zp->Fill(ZDChit->fLightPMQ);
 198                if(ZDChit->fVolume[1]==4)hPMQ4zp->Fill(ZDChit->fLightPMQ);
 199                EtotZP   += ZDChit->fEnergy;
 200                LtotZP   += (ZDChit->fLightPMQ) + (ZDChit->fLightPMC);
 201                LightCzp += ZDChit->fLightPMC;
 202                hspotzp->Fill(ZDChit->fXImpact,ZDChit->fYImpact);
 203              }
 204              if(ZDChit->fVolume[0]==3) { //ZEM
 205                hEzem->Fill(ZDChit->fEnergy);
 206                hPMzem->Fill(ZDChit->fLightPMC);
 207                hspotzem->Fill(ZDChit->fXImpact,ZDChit->fYImpact);
 208              }
 209
 210            }
 211            if(nhits!=0){
 212              hEzn->Fill(EtotZN);
 213              hLzn->Fill(LtotZN);
 214              hPMCzn->Fill(LightCzn);
 215              hEzp->Fill(EtotZP);
 216              hLzp->Fill(LtotZP);
 217              hPMCzp->Fill(LightCzp);
 218 //             printf("\n Histos var -> Ezn = %f, Lzn = %f, Ezp = %f, Lzp = %f \n\n",
 219 //                     EtotZN, LtotZN, EtotZP, LtotZP);
 220            }
 221          }//ZDC
 222       }//Track loop
 223    }//Hit loop
 224
 225 // Control prints
 226
 227
 228 if(detector == 1){ // ZN histos
 229    TCanvas *c1 = new TCanvas("c1","ZN hits",0,10,580,700);
 230    c1->cd();
 231    TPad *pad1= new TPad("pad1"," ",0.01,0.51,0.99,0.99);
 232    TPad *pad2= new TPad("pad2"," ",0.01,0.01,0.49,0.49);
 233    TPad *pad3= new TPad("pad3"," ",0.51,0.01,0.99,0.49);
 234    pad1->SetFillColor(18);
 235    pad2->SetFillColor(18);
 236    pad3->SetFillColor(18);
 237    pad1->Draw();
 238    pad2->Draw();
 239    pad3->Draw();
 240    pad1->cd();
 241    hspotzn->Draw();
 242    pad2->cd();
 243    hEzn->Draw();
 244    pad3->cd();
 245    hPMCzn->Draw();
 246
 247    TCanvas *c2 = new TCanvas("c2","ZN hits",600,10,600,700);
 248    c2->cd();
 249    TPad *pad4= new TPad("pad4"," ",0.01,0.51,0.49,0.99);
 250    TPad *pad5= new TPad("pad5"," ",0.51,0.51,0.99,0.99);
 251    TPad *pad6= new TPad("pad6"," ",0.01,0.01,0.49,0.49);
 252    TPad *pad7= new TPad("pad7"," ",0.51,0.01,0.99,0.49);
 253    pad4->SetFillColor(18);
 254    pad5->SetFillColor(18);
 255    pad6->SetFillColor(18);
 256    pad7->SetFillColor(18);
 257    pad4->Draw();
 258    pad5->Draw();
 259    pad6->Draw();
 260    pad7->Draw();
 261    pad4->cd();
 262    hPMQ1zn->Draw();
 263    pad5->cd();
 264    hPMQ2zn->Draw();
 265    pad6->cd();
 266    hPMQ3zn->Draw();
 267    pad7->cd();
 268    hPMQ4zn->Draw();
 269
 270    TCanvas *c3 = new TCanvas("c3","ZN hits",300,10,600,700);
 271    c3->cd();
 272    hLzn->Draw();
 273 }
 274
 275 if(detector == 2){ // ZP histos
 276    TCanvas *c1 = new TCanvas("c1","ZP hits",0,10,580,700);
 277    c1->cd();
 278    TPad *pad1= new TPad("pad1"," ",0.01,0.51,0.99,0.99);
 279    TPad *pad2= new TPad("pad2"," ",0.01,0.01,0.49,0.49);
 280    TPad *pad3= new TPad("pad3"," ",0.51,0.01,0.99,0.49);
 281    pad1->SetFillColor(18);
 282    pad2->SetFillColor(18);
 283    pad3->SetFillColor(18);
 284    pad1->Draw();
 285    pad2->Draw();
 286    pad3->Draw();
 287    pad1->cd();
 288    hspotzp->Draw();
 289    pad2->cd();
 290    hEzp->Draw();
 291    pad3->cd();
 292    hPMCzp->Draw();
 293
 294    TCanvas *c2 = new TCanvas("c2","ZP hits",600,10,600,700);
 295    c2->cd();
 296    TPad *pad4= new TPad("pad4"," ",0.01,0.51,0.49,0.99);
 297    TPad *pad5= new TPad("pad5"," ",0.51,0.51,0.99,0.99);
 298    TPad *pad6= new TPad("pad6"," ",0.01,0.01,0.49,0.49);
 299    TPad *pad7= new TPad("pad7"," ",0.51,0.01,0.99,0.49);
 300    pad4->SetFillColor(18);
 301    pad5->SetFillColor(18);
 302    pad6->SetFillColor(18);
 303    pad7->SetFillColor(18);
 304    pad4->Draw();
 305    pad5->Draw();
 306    pad6->Draw();
 307    pad7->Draw();
 308    pad4->cd();
 309    hPMQ1zp->Draw();
 310    pad5->cd();
 311    hPMQ2zp->Draw();
 312    pad6->cd();
 313    hPMQ3zp->Draw();
 314    pad7->cd();
 315    hPMQ4zp->Draw();
 316
 317    TCanvas *c3 = new TCanvas("c3","ZP hits",300,10,600,700);
 318    c3->cd();
 319    hLzp->Draw();
 320 }
 321
 322 if(detector == 3){ // ZEM histos
 323    TCanvas *c1 = new TCanvas("c1","ZEM hits",0,10,580,700);
 324    c1->cd();
 325    TPad *pad1= new TPad("pad1"," ",0.01,0.51,0.99,0.99);
 326    TPad *pad2= new TPad("pad2"," ",0.01,0.01,0.99,0.49);
 327    pad1->SetFillColor(18);
 328    pad2->SetFillColor(18);
 329    pad1->Draw();
 330    pad2->Draw();
 331    pad1->cd();
 332    hspotzem->Draw();
 333    pad2->cd();
 334    hEzem->Draw();
 335
 336    TCanvas *c2 = new TCanvas("c2","ZEM hits",600,10,600,700);
 337    c2->cd();
 338    hPMzem->Draw();
 339 }
 340 //   file->Close();
 341 }