ITS/ShowSPDRecPoints.C

   1 //////////////////////////////////////////////////////////////////
   2 // Macro to check clusters in the 2 SPD layers                  //
   3 // Provides:                                                    //
   4 //  2 canvases with                                             //
   5 //     - cluster loc and glob coordinates  (each layer)         //
   6 //  1 canvas with                                               //
   7 //      - cluster glob coordinates 2D and 3D                    //
   8 //  1 canvas with                                               //
   9 //      - correlations of #clusters for sectors                 //
  10 //      - correlations of #clusters for half-sectors            //
  11 //                                                              //
  12 //  Maria.Nicassio@ba.infn.it                                   //
  13 //  Domenico.Elia@ba.infn.it                                    //
  14 //////////////////////////////////////////////////////////////////
  15
  16 #if !defined(__CINT__) || defined(__MAKECINT__)
  17
  18 #include <Riostream.h>
  19 #include <TFile.h>
  20 #include <TTree.h>
  21 #include <TBranch.h>
  22 #include <TCanvas.h>
  23 #include <TClonesArray.h>
  24 #include <TH1.h>
  25 #include <TH2.h>
  26 #include <TH3.h>
  27 #include <TStyle.h>
  28 #include <TProfile.h>
  29 #include <TClassTable.h>
  30 #include <TInterpreter.h>
  31 #include <TGraph2D.h>
  32 #include <TGraph.h>
  33 #include <TGeoManager.h>
  34
  35 #include "AliCDBManager.h"
  36 #include "AliRunLoader.h"
  37 #include "AliESD.h"
  38 #include "AliRun.h"
  39 #include "AliGeomManager.h"
  40 #include "AliITS.h"
  41 #include "AliITSgeom.h"
  42 #include "AliITSLoader.h"
  43 #include <AliITSRecPoint.h>
  44 #include <AliHeader.h>
  45 #include <AliITSDetTypeRec.h>
  46
  47 #endif
  48
  49 /* $Id$ */
  50
  51 void ShowSPDRecPoints(Int_t RunStart, Int_t RunStop){
  52
  53   // Set data directory
  54   Char_t str[256];
  55   Char_t* dir = "/data/alipix/PhysicsEvents/pp/ppProd/pdc07/oldgeom";
  56
  57   // Variables for histo booking and filling
  58   Int_t modmin=0;
  59   Int_t modmax=240;
  60   Int_t totmod=modmax-modmin;
  61
  62   Float_t xlim[2]={4.5,8.};
  63   Float_t zlim[2]={15.,15.};
  64
  65   Int_t nClusters[2];
  66
  67   Float_t clustersXCoord[2][200];
  68   Float_t clustersYCoord[2][200];
  69   Float_t clustersZCoord[2][200];
  70
  71   Float_t cluGlo[3]={0.,0.,0.};
  72
  73   Int_t nClustersPerLayerPerSector[2][10];
  74   Int_t nClustersPerLayerPerHalfSector[2][20];
  75
  76   Int_t iSector=0;
  77   Int_t iHalfSector=0;
  78
  79   // Booking of histograms
  80   gStyle->SetPalette(1,0);
  81   TH1F* hlayer= new TH1F("hlayer","",6,0.,6.);
  82   TH1F** hmod = new TH1F*[2];
  83   TH1F** hxl  = new TH1F*[2];
  84   TH1F** hzl  = new TH1F*[2];
  85   TH1F** hxg  = new TH1F*[2];
  86   TH1F** hyg  = new TH1F*[2];
  87   TH1F** hzg  = new TH1F*[2];
  88   TH1F** hr   = new TH1F*[2];
  89   TH1F** hphi = new TH1F*[2];
  90   TH1F** hMultSPDcl = new TH1F*[2];
  91   TH1F** hType = new TH1F*[2];  // cluster type ?
  92
  93   TH2F** hr_phi   = new TH2F*[2];
  94   TH2F** hx_y   = new TH2F*[2];
  95   TH3F** hx_y_z   = new TH3F*[2];
  96
  97   TH2F** hnCl2_nCl1_Sec = new TH2F*[10];
  98   TProfile** hnCl2vsnCl1_Sec = new TProfile*[10];
  99   TH2F** hnCl2_nCl1_HSec = new TH2F*[20];
 100 //  TProfile** hnCl2vsnCl1_HSec = new TProfile*[20];
 101   TH2F** hNy_Nz = new TH2F*[2];  // y and z length
 102   TH2F** hPhi_Z = new TH2F*[2];
 103
 104   TH2F *hMultSPDcl2_MultSPDcl1 =
 105             new TH2F("nCl2_nCl1","# SPD clusters on layer 2 vs # on layer 1",200,0.,200.,200,0.,200.);
 106   hMultSPDcl2_MultSPDcl1->GetXaxis()->SetTitle("# clusters (inner layer)");
 107   hMultSPDcl2_MultSPDcl1->GetYaxis()->SetTitle("# clusters (outer layer)");
 108   Char_t name[10];
 109   Char_t title[20];
 110   for (Int_t iLay=0;iLay<2;iLay++) {
 111     sprintf(name,"hmod%d",iLay+1);
 112     hmod[iLay]=new TH1F(name,"SPD clusters - Module number",totmod,modmin,modmax);
 113     hmod[iLay]->GetXaxis()->SetTitle("Module number");
 114     hmod[iLay]->GetYaxis()->SetTitle("Entries");
 115     sprintf(name,"hxloc%d",iLay+1);
 116     hxl[iLay]=new TH1F(name,"SPD clusters - Local x coordinate",100,-4.,4.);
 117     hxl[iLay]->GetXaxis()->SetTitle("Local x [cm]");
 118     hxl[iLay]->GetYaxis()->SetTitle("Entries");
 119     sprintf(name,"hzloc%d",iLay+1);
 120     hzl[iLay]=new TH1F(name,"SPD clusters - Local z coordinate",100,-4.,4.);
 121     hzl[iLay]->GetXaxis()->SetTitle("Local z [cm]");
 122     hzl[iLay]->GetYaxis()->SetTitle("Entries");
 123     sprintf(name,"hxgl%d",iLay+1);
 124     hxg[iLay]=new TH1F(name,"SPD clusters - Global x coordinate",100,-xlim[iLay],xlim[iLay]);
 125     hxg[iLay]->GetXaxis()->SetTitle("Global x [cm]");
 126     hxg[iLay]->GetYaxis()->SetTitle("Entries");
 127     sprintf(name,"hygl%d",iLay+1);
 128     hyg[iLay]=new TH1F(name,"SPD clusters - Global y coordinate",100,-xlim[iLay],xlim[iLay]);
 129     hyg[iLay]->GetXaxis()->SetTitle("Global y [cm]");
 130     hyg[iLay]->GetYaxis()->SetTitle("Entries");
 131     sprintf(name,"hzgl%d",iLay+1);
 132     hzg[iLay]=new TH1F(name,"SPD clusters - Global z coordinate",150,-zlim[iLay],zlim[iLay]);
 133     hzg[iLay]->GetXaxis()->SetTitle("Global z [cm]");
 134     hzg[iLay]->GetYaxis()->SetTitle("Entries");
 135     sprintf(name,"hr%d",iLay+1);
 136     hr[iLay]=new TH1F(name,"SPD clusters - r",100,0.,50.);
 137     hr[iLay]->GetXaxis()->SetTitle("r [cm]");
 138     hr[iLay]->GetYaxis()->SetTitle("Entries");
 139     sprintf(name,"hphi%d",iLay+1);
 140     hphi[iLay]=new TH1F(name,"SPD clusters - #phi",100,0.,2*TMath::Pi());
 141     hphi[iLay]->GetXaxis()->SetTitle("#phi [rad]");
 142     hphi[iLay]->GetYaxis()->SetTitle("Entries");
 143     sprintf(name,"hType%d",iLay+1);
 144     hType[iLay]=new TH1F(name,"SPD clusters - Type",100,0.,100.);
 145     hType[iLay]->GetXaxis()->SetTitle("Cluster type");
 146     hType[iLay]->GetYaxis()->SetTitle("Entries");
 147
 148     sprintf(name,"hMultSPDcl%d",iLay+1);
 149     hMultSPDcl[iLay]=new TH1F(name,"Cluster multiplicity",200,0.,200.);
 150     hMultSPDcl[iLay]->GetXaxis()->SetTitle("Cluster multiplicity");
 151     hMultSPDcl[iLay]->GetYaxis()->SetTitle("Entries");
 152
 153     sprintf(name,"hNy_Nz%d",iLay+1);
 154     hNy_Nz[iLay]=new TH2F(name,"SPD clusters - Length",100,0.,100.,100,0.,100.);
 155     hNy_Nz[iLay]->GetXaxis()->SetTitle("z length");
 156     hNy_Nz[iLay]->GetYaxis()->SetTitle("y length");
 157     sprintf(name,"hPhi_Z%d",iLay+1);
 158     hPhi_Z[iLay]=new TH2F(name,"SPD clusters - #phi vs z",150,-zlim[iLay],zlim[iLay],100,0.,2*TMath::Pi());
 159     hPhi_Z[iLay]->GetXaxis()->SetTitle("Global z [cm]");
 160     hPhi_Z[iLay]->GetYaxis()->SetTitle("#phi [rad]");
 161     sprintf(name,"hr_phi%d",iLay+1);
 162     hr_phi[iLay]=new TH2F(name,"SPD clusters - #phi_r",500,0.,50.,100,0.,2*TMath::Pi());
 163     hr_phi[iLay]->GetXaxis()->SetTitle("r [cm]");
 164     hr_phi[iLay]->GetYaxis()->SetTitle("#phi [rad]");
 165     sprintf(name,"hx_y%d",iLay+1);
 166     hx_y[iLay]=new TH2F(name,"SPD clusters - y_x",200,-10.,10.,200,-10.,10.);
 167     hx_y[iLay]->GetXaxis()->SetTitle("x [cm]");
 168     hx_y[iLay]->GetYaxis()->SetTitle("y [cm]");
 169     sprintf(name,"hx_y_z%d",iLay+1);
 170     hx_y_z[iLay]=new TH3F(name,"SPD clusters - y_x",200,-10.,10.,200,-10.,10.,150,-15.,15.);
 171     hx_y_z[iLay]->GetXaxis()->SetTitle("z [cm]");
 172     hx_y_z[iLay]->GetYaxis()->SetTitle("x [cm]");
 173     hx_y_z[iLay]->GetZaxis()->SetTitle("y [cm]");
 174   }
 175   for (Int_t iSec=0; iSec<10; iSec++) {
 176     sprintf(name,"hnCl2_nCl1_Sector%d",iSec);
 177     sprintf(title,"Sector %d",iSec+1);
 178     hnCl2_nCl1_Sec[iSec]=new TH2F(name,title,200,0.,200.,200,0.,200.);
 179     hnCl2_nCl1_Sec[iSec]->GetXaxis()->SetTitle("# clusters (inner layer)");
 180     hnCl2_nCl1_Sec[iSec]->GetYaxis()->SetTitle("# clusters (outer layer)");
 181     sprintf(name,"hnCl2vsnCl1_Sector%d",iSec);
 182     sprintf(title,"Sector %d",iSec+1);
 183     hnCl2vsnCl1_Sec[iSec]=new TProfile(name,title,200,0.,200.,0.,200.);
 184     hnCl2vsnCl1_Sec[iSec]->GetXaxis()->SetTitle("# clusters (inner layer)");
 185     hnCl2vsnCl1_Sec[iSec]->GetYaxis()->SetTitle("# clusters (outer layer)");
 186   }
 187   for (Int_t iHalfSec=0; iHalfSec<20; iHalfSec++) {
 188     sprintf(name,"hnCl2_nCl1_HalfSector%d",iHalfSec);
 189     sprintf(title,"Half-Sector %d",iHalfSec+1);
 190     hnCl2_nCl1_HSec[iHalfSec]=new TH2F(name,title,200,0.,200.,200,0.,200.);
 191     hnCl2_nCl1_HSec[iHalfSec]->GetXaxis()->SetTitle("# clusters (inner layer)");
 192     hnCl2_nCl1_HSec[iHalfSec]->GetYaxis()->SetTitle("# clusters (outer layer)");
 193   }
 194
 195   // Loop over runs
 196   for (Int_t run=RunStart; run<RunStop+1; run++) {
 197
 198     // setup galice and runloader
 199 //    cout << "File nr --> " << run << endl;
 200
 201     if (gClassTable->GetID("AliRun") < 0) {
 202       gInterpreter->ExecuteMacro("loadlibs.C");
 203     }
 204     else {
 205       if (gAlice){
 206         delete AliRunLoader::Instance();
 207         delete gAlice;
 208         gAlice=0;
 209       }
 210     }
 211
 212     // Set OfflineConditionsDataBase if needed
 213     AliCDBManager* man = AliCDBManager::Instance();
 214     if (!man->IsDefaultStorageSet()) {
 215       printf("Setting a local default storage and run number 0\n");
 216       man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
 217       man->SetRun(0);
 218     }
 219     else {
 220       printf("Using deafult storage \n");
 221     }
 222
 223     // retrives geometry
 224     if (!gGeoManager) {
 225       sprintf(str,"%s/ppMinBias%04d/geometry.root",dir,run);
 226       AliGeomManager::LoadGeometry(str);
 227     }
 228
 229     sprintf(str,"%s/ppMinBias%04d/galice.root",dir,run);
 230     AliRunLoader*  rl = AliRunLoader::Open(str);
 231
 232     if (rl == 0x0){
 233       cerr<<"Can not open session RL=NULL"<< endl;
 234       return;
 235     }
 236     Int_t retval = rl->LoadgAlice();
 237     if (retval){
 238       cerr<<"LoadgAlice returned error"<<endl;
 239       return;
 240     }
 241     gAlice=rl->GetAliRun();
 242
 243     retval = rl->LoadHeader();
 244     if (retval){
 245       cerr<<"LoadHeader returned error"<<endl;
 246       return;
 247     }
 248
 249     AliITSLoader* ITSloader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
 250     if (!ITSloader){
 251       cerr<<"ITS loader not found"<<endl;
 252       return;
 253     }
 254     ITSloader->LoadRecPoints("read");
 255
 256     AliITS *ITS  = (AliITS*)gAlice->GetModule("ITS");
 257     ITS->SetTreeAddress();
 258     AliITSgeom *geom = ITS->GetITSgeom();
 259     if (!geom) {
 260       cout << " Can't get the ITS geometry!" << endl;
 261       return ;
 262     }
 263
 264     AliITSDetTypeRec* detTypeRec = new AliITSDetTypeRec();
 265
 266     detTypeRec->SetITSgeom(ITSloader->GetITSgeom());
 267     detTypeRec->SetDefaults();
 268
 269     Int_t nEvents=rl->GetNumberOfEvents();
 270     printf("Total Number of events = %d\n",nEvents);
 271
 272     for (Int_t iev=0;iev<nEvents;iev++){
 273       rl->GetEvent(iev);
 274
 275       TTree *TR = ITSloader->TreeR();
 276       TClonesArray* ITSClusters;
 277       ITSClusters  = detTypeRec->RecPoints();
 278       TBranch *branch = 0;
 279       if (TR && ITSClusters) {
 280         branch = ITSloader->TreeR()->GetBranch("ITSRecPoints");
 281         if (branch) branch->SetAddress(&ITSClusters);
 282       }
 283
 284 //      Int_t nparticles = rl->GetHeader()->GetNtrack();
 285 //      cout<<"Event # "<<iev<<"   #Particles="<<nparticles<<endl;
 286
 287       // Reset cluster counters
 288       for (Int_t iLay=0; iLay<2; iLay++) {
 289         nClusters[iLay]=0;
 290         for (Int_t iSec=0; iSec<10; iSec++) {
 291           nClustersPerLayerPerSector[iLay][iSec]=0;
 292         }
 293         for (Int_t iHSec=0; iHSec<20; iHSec++) {
 294           nClustersPerLayerPerHalfSector[iLay][iHSec]=0;
 295         }
 296       }
 297       for (Int_t subd=0;subd<1;subd++) {
 298
 299         Int_t first = geom->GetStartDet(subd);
 300         Int_t last = geom->GetLastDet(subd);
 301         for (Int_t mod=first; mod<=last; mod++) {
 302           detTypeRec->ResetRecPoints();
 303           branch->GetEvent(mod);
 304           Int_t nrecp = ITSClusters->GetEntries();
 305           while (nrecp--) {
 306
 307             AliITSRecPoint *recp = (AliITSRecPoint*)ITSClusters->At(nrecp);
 308
 309             Int_t lay=recp->GetLayer();
 310 //            cout<<"lay"<<lay<<endl;
 311             recp->GetGlobalXYZ(cluGlo);
 312
 313             Float_t rad=TMath::Sqrt(cluGlo[0]*cluGlo[0]+cluGlo[1]*cluGlo[1]);
 314             Float_t phi= TMath::Pi() + TMath::ATan2(-cluGlo[1],-cluGlo[0]);
 315             hlayer->Fill(lay);
 316             hmod[lay]->Fill(mod);
 317             hzl[lay]->Fill(recp->GetDetLocalZ());
 318             hxl[lay]->Fill(recp->GetDetLocalX());
 319             hzg[lay]->Fill(cluGlo[2]);
 320             hyg[lay]->Fill(cluGlo[1]);
 321             hxg[lay]->Fill(cluGlo[0]);
 322             hr[lay]->Fill(rad);
 323             hphi[lay]->Fill(phi);
 324             hPhi_Z[lay]->Fill(cluGlo[2],phi);
 325             hr_phi[lay]->Fill(rad,phi);
 326             hx_y[lay]->Fill(cluGlo[0],cluGlo[1]);
 327             hx_y_z[lay]->Fill(cluGlo[2],cluGlo[0],cluGlo[1]);
 328             clustersXCoord[lay][nClusters[lay]]=cluGlo[0];
 329             clustersYCoord[lay][nClusters[lay]]=cluGlo[1];
 330             clustersZCoord[lay][nClusters[lay]]=cluGlo[2];
 331             nClusters[lay]++;
 332
 333             hType[lay]->Fill(recp->GetType());
 334 //            cout<<"Clusters type"<<recp->GetType()<<endl;
 335             hNy_Nz[lay]->Fill(recp->GetNz(),recp->GetNy());
 336
 337             // Set Sector number and increase the counter
 338             if (lay==0) {
 339               for (Int_t nRange=0; nRange<10; nRange++) {
 340                 if ((mod>=nRange*8) && (mod<=(nRange*8+7))) {
 341                   iSector = nRange;
 342                   nClustersPerLayerPerSector[lay][iSector]++;
 343                   break;
 344                 }
 345               }
 346             }
 347             if (lay==1) {
 348               for (Int_t nRange=0; nRange<10; nRange++) {
 349                 if ((mod>=80+nRange*16) && (mod<=(80+nRange*16+15))) {
 350                   iSector = nRange;
 351                   nClustersPerLayerPerSector[lay][iSector]++;
 352                   break;
 353                 }
 354               }
 355             }
 356             // Set HalfSector number and increase the counter
 357             if (lay==0) {
 358               for (Int_t nRange=0; nRange<20; nRange++) {
 359                  if ((mod>=nRange*4) && (mod<=(nRange*4+3))) {
 360                    iHalfSector = nRange;
 361                    nClustersPerLayerPerHalfSector[lay][iHalfSector]++;
 362                    break;
 363                  }
 364               }
 365             } else {
 366               for (Int_t nRange=0; nRange<20; nRange++) {
 367                 if ((mod>=80+nRange*4) && (mod<=(80+nRange*4+3))) {
 368                   iHalfSector = nRange;
 369                   nClustersPerLayerPerHalfSector[lay][iHalfSector]++;
 370                   break;
 371                 }
 372               }
 373             }
 374           }
 375         }
 376       }
 377
 378       for (Int_t iLay=0; iLay<2; iLay++)   hMultSPDcl[iLay]->Fill(nClusters[iLay]);
 379       for (Int_t iSec=0; iSec<10; iSec++) {
 380         hnCl2_nCl1_Sec[iSec]->Fill(nClustersPerLayerPerSector[0][iSec],nClustersPerLayerPerSector[1][iSec]);
 381         hnCl2vsnCl1_Sec[iSec]->Fill(nClustersPerLayerPerSector[0][iSec],nClustersPerLayerPerSector[1][iSec]);
 382       }
 383       for (Int_t iHSec=0; iHSec<20; iHSec++) {
 384         hnCl2_nCl1_HSec[iHSec]->Fill(nClustersPerLayerPerHalfSector[0][iHSec],nClustersPerLayerPerHalfSector[1][iHSec]);
 385       }
 386       hMultSPDcl2_MultSPDcl1->Fill(nClusters[0],nClusters[1]);
 387
 388     } //end loop over events
 389     rl->UnloadAll();
 390     delete rl;
 391
 392   } //end loop over runs
 393
 394   // Draw and Write histos
 395
 396   TFile* fout = new TFile("out_ShowSPDRecPoints.root","RECREATE");
 397
 398   hlayer->Write();
 399 //  cev0->Write();
 400
 401   TCanvas **c=new TCanvas*[2];
 402   Char_t ctit[30];
 403   for(Int_t iLay=0;iLay<2;iLay++){
 404     hNy_Nz[iLay]->Write();
 405     sprintf(name,"can%d",iLay+1);
 406     sprintf(ctit,"Layer %d",iLay+1);
 407     c[iLay]=new TCanvas(name,ctit,1200,900);
 408     c[iLay]->Divide(3,3);
 409     c[iLay]->cd(1);
 410     hmod[iLay]->Draw();
 411     hmod[iLay]->Write();
 412     c[iLay]->cd(2);
 413     hxl[iLay]->Draw();
 414     hxl[iLay]->Write();
 415     c[iLay]->cd(3);
 416     hzl[iLay]->Draw();
 417     hzl[iLay]->Write();
 418     c[iLay]->cd(4);
 419     hxg[iLay]->Draw();
 420     hxg[iLay]->Write();
 421     c[iLay]->cd(5);
 422     hyg[iLay]->Draw();
 423     hyg[iLay]->Write();
 424     c[iLay]->cd(6);
 425     hzg[iLay]->Draw();
 426     hzg[iLay]->Write();
 427     c[iLay]->cd(7);
 428     hr[iLay]->Draw();
 429     hr[iLay]->Write();
 430     c[iLay]->cd(8);
 431     hphi[iLay]->Draw();
 432     hphi[iLay]->Write();
 433     c[iLay]->cd(9);
 434     hPhi_Z[iLay]->Draw("colz");
 435     hPhi_Z[iLay]->Write();
 436     hr_phi[iLay]->Write();
 437     hx_y[iLay]->Write();
 438     hx_y_z[iLay]->Write();
 439   }
 440
 441   TCanvas *cCoord=new TCanvas("cCoord","Cluster coordinates",1200,900);
 442   cCoord->Divide(2,2);
 443
 444   for (Int_t iLay=0;iLay<2;iLay++) {
 445     cCoord->cd(1);
 446     hx_y[iLay]->SetMarkerStyle(22);
 447     hx_y[iLay]->SetMarkerSize(0.3);
 448     hx_y[iLay]->SetMarkerColor(iLay+1);
 449     if (iLay==0) hx_y[iLay]->Draw("p");
 450     else hx_y[iLay]->Draw("p,same");
 451     cCoord->cd(2);
 452     hx_y_z[iLay]->SetMarkerStyle(23);
 453     hx_y_z[iLay]->SetMarkerSize(0.3);
 454     hx_y_z[iLay]->SetMarkerColor(iLay+1);
 455     hx_y_z[iLay]->Draw("p,same");
 456     cCoord->cd(3);
 457     hr_phi[iLay]->SetMarkerColor(iLay+1);
 458     hr_phi[iLay]->Draw("p,same");
 459   }
 460   TCanvas *cCorrelations_Sectors=new TCanvas("cCorrelations_S","SPD cluster correlations (Sectors)",1200,900);
 461   cCorrelations_Sectors->Divide(3,5);
 462   cCorrelations_Sectors->cd(1);
 463   hMultSPDcl2_MultSPDcl1->Draw();
 464   hMultSPDcl2_MultSPDcl1->Write();
 465   cCorrelations_Sectors->cd(2);
 466   hMultSPDcl[0]->Draw();
 467   hMultSPDcl[0]->Write();
 468   cCorrelations_Sectors->cd(3);
 469   hMultSPDcl[1]->Draw();
 470   hMultSPDcl[1]->Write();
 471   for (Int_t iS=0; iS<10; ++iS) {
 472     cCorrelations_Sectors->cd(iS+4);
 473     hnCl2_nCl1_Sec[iS]->Draw();
 474     hnCl2_nCl1_Sec[iS]->Write();
 475     hnCl2vsnCl1_Sec[iS]->Write();
 476   }
 477
 478   TCanvas *cCorrelations_HSectors=new TCanvas("cCorrelations_HS","SPD cluster correlations (Half-Sectors)",1200,900);
 479   cCorrelations_HSectors->Divide(5,4);
 480
 481   for (Int_t iHS=0; iHS<20; ++iHS) {
 482     cCorrelations_HSectors->cd(iHS+1);
 483     hnCl2_nCl1_HSec[iHS]->Write();
 484     hnCl2_nCl1_HSec[iHS]->Draw();
 485   }
 486
 487   fout->Close();
 488
 489 }