MFT/MergeClustersMFT.C

   1 #include "AliMFTCluster.h"
   2 #include "TFile.h"
   3 #include "AliMFTCluster.h"
   4 #include "TObjArray.h"
   5 #include "TClonesArray.h"
   6 #include "TTree.h"
   7 #include "AliMFTConstants.h"
   8 #include "TBranch.h"
   9 #include "TRandom.h"
  10 #include "TMath.h"
  11
  12 TFile *fileInPileUp=0, *fileInSignal, *fileInUnderlying, *fileOutPileUp=0, *fileOutTotal;
  13
  14 enum {kMergePileUpClusters_pp, kMergePileUpClusters_PbPb, kMergeAllClusters_pp, kMergeAllClusters_PbPb};
  15 enum {k_pp, k_PbPb};
  16
  17 const Char_t *fPileUpDir[2] = {"$UnderlyingDirPP/spreadVtx",      // local directory with pile-up event p-p
  18                                "$UnderlyingDirPbPb/spreadVtx"};   // local directory with pile-up event Pb-Pb
  19
  20 const Int_t fNRunsPileUp_pp   = 100;
  21 const Int_t fNRunsPileUp_PbPb =   3;
  22
  23 const Int_t fNRunsPileUpAvailable[2] = {500, 500};    // pp, PbPb
  24
  25 TClonesArray *fRecPointsPerPlaneIn[AliMFTConstants::fNMaxPlanes] = {0};
  26 TClonesArray *fRecPointsPerPlaneOut[AliMFTConstants::fNMaxPlanes] = {0};
  27
  28 void MergePileUpClusters(Int_t *runs, Int_t nRunsPileUp, const Char_t *pileUpDir, Int_t nEvents);
  29 void MergeAllClusters();
  30
  31 //====================================================================================================================================================
  32
  33 void MergeClustersMFT(Int_t option = 9999,
  34                       Int_t nEvents = 9999,
  35                       Int_t maxPileUp = 9999) {
  36
  37   if (option==kMergePileUpClusters_pp) {
  38
  39     Int_t runs_pp[fNRunsPileUp_pp] = {0};
  40     Int_t nFiles_pp = 0;
  41     while (nFiles_pp<TMath::Min(fNRunsPileUp_pp,maxPileUp)) {
  42       Int_t run = gRandom->Integer(fNRunsPileUpAvailable[k_pp]);
  43       fileInPileUp = new TFile(Form("%s/run_%d/MFT.RecPoints.MCShifted.root",fPileUpDir[k_pp],run));
  44       if (!fileInPileUp->IsOpen() || !fileInPileUp->cd("Event0")) {
  45         delete fileInPileUp;
  46         continue;
  47       }
  48       Bool_t runAlreadyLoaded=kFALSE;
  49       for (Int_t iAddedRun=0; iAddedRun<nFiles_pp; iAddedRun++) {
  50         if (run==runs_pp[iAddedRun]) {
  51           runAlreadyLoaded=kTRUE;
  52           break;
  53         }
  54       }
  55       if (!runAlreadyLoaded) runs_pp[nFiles_pp++] = run;
  56       delete fileInPileUp;
  57     }
  58     MergePileUpClusters(runs_pp, TMath::Min(fNRunsPileUp_pp,maxPileUp), fPileUpDir[k_pp], nEvents);
  59   }
  60
  61   else if (option==kMergePileUpClusters_PbPb) {
  62     Int_t runs_PbPb[fNRunsPileUp_PbPb] = {0};
  63     Int_t nFiles_PbPb = 0;
  64     while (nFiles_PbPb<TMath::Min(fNRunsPileUp_PbPb,maxPileUp)) {
  65       Int_t run = gRandom->Integer(fNRunsPileUpAvailable[k_PbPb]);
  66       fileInPileUp = new TFile(Form("%s/run_%d/MFT.RecPoints.MCShifted.root",fPileUpDir[k_PbPb],run));
  67       if (!fileInPileUp->IsOpen() || !fileInPileUp->cd("Event0")) {
  68         delete fileInPileUp;
  69         continue;
  70       }
  71       Bool_t runAlreadyLoaded=kFALSE;
  72       for (Int_t iAddedRun=0; iAddedRun<nFiles_PbPb; iAddedRun++) {
  73         if (run==runs_PbPb[iAddedRun]) {
  74           runAlreadyLoaded=kTRUE;
  75           break;
  76         }
  77       }
  78       if (!runAlreadyLoaded) runs_PbPb[nFiles_PbPb++] = run;
  79       delete fileInPileUp;
  80     }
  81     MergePileUpClusters(runs_PbPb, TMath::Min(fNRunsPileUp_PbPb,maxPileUp), fPileUpDir[k_PbPb], nEvents);
  82   }
  83
  84   if (option>=kMergeAllClusters_pp) MergeAllClusters();
  85
  86 }
  87
  88 //====================================================================================================================================================
  89
  90 void MergePileUpClusters(Int_t *runs, Int_t nRunsPileUp, const Char_t *pileUpDir, Int_t nEvents) {
  91
  92   fileOutPileUp = new TFile("./MFT.RecPoints.PileUp.root", "recreate");
  93
  94   Bool_t eventExist = kTRUE;
  95   Int_t iEv=0;
  96
  97   while (eventExist && iEv<nEvents) {
  98
  99     for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 100       fRecPointsPerPlaneOut[iPlane] = new TClonesArray("AliMFTCluster");
 101     }
 102
 103     if (iEv) fileOutPileUp = new TFile("./MFT.RecPoints.PileUp.root", "update");
 104
 105     printf("Merging Event %d\n",iEv);
 106
 107     fileOutPileUp-> cd();
 108     TTree *treeOutPileUp = new TTree("TreeR_OutPileUp", "Reconstructed Points Container");
 109
 110     Int_t nFilesAdded = 0;
 111     while (nFilesAdded<nRunsPileUp) {
 112       fileInPileUp = new TFile(Form("%s/run_%d/MFT.RecPoints.MCShifted.root",pileUpDir,runs[nFilesAdded++]));
 113       printf("Merging file %s\n", fileInPileUp->GetName());
 114       if (!fileInPileUp->cd(Form("Event%d",iEv))) {
 115         eventExist = kFALSE;
 116         break;
 117       }
 118       for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 119         fRecPointsPerPlaneIn[iPlane] = new TClonesArray("AliMFTCluster");
 120       }
 121       TTree *treeInPileUp = (TTree*) gDirectory->Get("TreeR");
 122       treeInPileUp -> SetName("TreeR_InPileUp");
 123       Int_t iPlane=0;
 124       while (treeInPileUp->GetBranch(Form("Plane_%02d",iPlane))) {
 125         //      printf("Plane %02d\n",iPlane);
 126         treeInPileUp ->SetBranchAddress(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneIn[iPlane]));
 127         TBranch *branch = treeOutPileUp->GetBranch(Form("Plane_%02d",iPlane));
 128         if (!branch) treeOutPileUp->Branch(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneOut[iPlane]));
 129         iPlane++;
 130       }
 131       iPlane=0;
 132       treeInPileUp -> GetEntry(0);
 133       while (treeInPileUp->GetBranch(Form("Plane_%02d",iPlane))) {
 134         Int_t nClusters = fRecPointsPerPlaneIn[iPlane]->GetEntries();
 135         for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
 136           //      printf("Cluster %4d\n",iCluster);
 137           AliMFTCluster *newCluster = (AliMFTCluster*) fRecPointsPerPlaneIn[iPlane]->At(iCluster);
 138           new ((*fRecPointsPerPlaneOut[iPlane])[fRecPointsPerPlaneOut[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
 139         }
 140         iPlane++;
 141       }
 142       for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 143         delete fRecPointsPerPlaneIn[jPlane];
 144       }
 145       delete treeInPileUp;
 146       fileInPileUp -> Close();
 147       delete fileInPileUp;
 148     }
 149     if (eventExist) {
 150       treeOutPileUp -> Fill();
 151       fileOutPileUp -> mkdir(Form("Event%d",iEv));
 152       fileOutPileUp -> cd(Form("Event%d",iEv));
 153       treeOutPileUp -> SetName("TreeR");
 154       treeOutPileUp -> Write();
 155       delete treeOutPileUp;
 156       fileOutPileUp -> Close();
 157       delete fileOutPileUp;
 158       for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 159         delete fRecPointsPerPlaneOut[iPlane];
 160       }
 161       iEv++;
 162     }
 163   }
 164
 165 }
 166
 167 //====================================================================================================================================================
 168
 169 void MergeAllClusters() {
 170
 171   fileOutTotal = new TFile("./MFT.RecPoints.root", "recreate");
 172
 173   Bool_t eventExist = kTRUE;
 174   Int_t iEv=0;
 175
 176   while (eventExist) {
 177
 178     for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 179       fRecPointsPerPlaneOut[iPlane] = new TClonesArray("AliMFTCluster");
 180     }
 181
 182     if (iEv) fileOutTotal = new TFile("./MFT.RecPoints.root", "update");
 183
 184     printf("Merging Event %d\n",iEv);
 185
 186     fileOutTotal-> cd();
 187     TTree *treeOutTotal = new TTree("TreeR_Out", "Reconstructed Points Container");
 188
 189     // Adding Pile-Up Clusters
 190
 191     fileInPileUp = new TFile("./MFT.RecPoints.PileUp.root");
 192     printf("Merging file %s\n", fileInPileUp->GetName());
 193     if (!fileInPileUp->cd(Form("Event%d",iEv))) {
 194       eventExist = kFALSE;
 195       break;
 196     }
 197     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 198       fRecPointsPerPlaneIn[jPlane] = new TClonesArray("AliMFTCluster");
 199     }
 200     TTree *treeInPileUp = (TTree*) gDirectory->Get("TreeR");
 201     treeInPileUp -> SetName("TreeR_InPileUp");
 202     Int_t iPlane=0;
 203     while (treeInPileUp->GetBranch(Form("Plane_%02d",iPlane))) {
 204       //        printf("Plane %02d\n",iPlane);
 205       treeInPileUp ->SetBranchAddress(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneIn[iPlane]));
 206       TBranch *branch = treeOutTotal->GetBranch(Form("Plane_%02d",iPlane));
 207       if (!branch) treeOutTotal->Branch(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneOut[iPlane]));
 208       iPlane++;
 209     }
 210     iPlane=0;
 211     treeInPileUp -> GetEntry(0);
 212     while (treeInPileUp->GetBranch(Form("Plane_%02d",iPlane))) {
 213       Int_t nClusters = fRecPointsPerPlaneIn[iPlane]->GetEntries();
 214       for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
 215         //        printf("Cluster %4d\n",iCluster);
 216         AliMFTCluster *newCluster = (AliMFTCluster*) fRecPointsPerPlaneIn[iPlane]->At(iCluster);
 217         new ((*fRecPointsPerPlaneOut[iPlane])[fRecPointsPerPlaneOut[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
 218       }
 219       iPlane++;
 220     }
 221     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 222       delete fRecPointsPerPlaneIn[jPlane];
 223     }
 224     delete treeInPileUp;
 225     fileInPileUp -> Close();
 226     delete fileInPileUp;
 227
 228     // Adding Underlying Event Clusters
 229
 230     fileInUnderlying = new TFile("./MFT.RecPoints.Underlying.root");
 231     printf("Merging file %s\n", fileInUnderlying->GetName());
 232     if (!fileInUnderlying->cd(Form("Event%d",iEv))) {
 233       eventExist = kFALSE;
 234       break;
 235     }
 236     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 237       fRecPointsPerPlaneIn[jPlane] = new TClonesArray("AliMFTCluster");
 238     }
 239     TTree *treeInUnderlying = (TTree*) gDirectory->Get("TreeR");
 240     treeInUnderlying -> SetName("TreeR_InUnderlying");
 241     iPlane=0;
 242     while (treeInUnderlying->GetBranch(Form("Plane_%02d",iPlane))) {
 243       //        printf("Plane %02d\n",iPlane);
 244       treeInUnderlying ->SetBranchAddress(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneIn[iPlane]));
 245       TBranch *branch = treeOutTotal->GetBranch(Form("Plane_%02d",iPlane));
 246       if (!branch) treeOutTotal->Branch(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneOut[iPlane]));
 247       iPlane++;
 248     }
 249     iPlane=0;
 250     treeInUnderlying -> GetEntry(0);
 251     while (treeInUnderlying->GetBranch(Form("Plane_%02d",iPlane))) {
 252       Int_t nClusters = fRecPointsPerPlaneIn[iPlane]->GetEntries();
 253       for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
 254         //        printf("Cluster %4d\n",iCluster);
 255         AliMFTCluster *newCluster = (AliMFTCluster*) fRecPointsPerPlaneIn[iPlane]->At(iCluster);
 256         new ((*fRecPointsPerPlaneOut[iPlane])[fRecPointsPerPlaneOut[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
 257       }
 258       iPlane++;
 259     }
 260     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 261       delete fRecPointsPerPlaneIn[jPlane];
 262     }
 263     delete treeInUnderlying;
 264     fileInUnderlying -> Close();
 265     delete fileInUnderlying;
 266
 267     // Adding Signal Clusters
 268
 269     fileInSignal = new TFile("./MFT.RecPoints.Signal.root");
 270     printf("Merging file %s\n", fileInSignal->GetName());
 271     if (!fileInSignal->cd(Form("Event%d",iEv))) {
 272       eventExist = kFALSE;
 273       break;
 274     }
 275     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 276       fRecPointsPerPlaneIn[jPlane] = new TClonesArray("AliMFTCluster");
 277     }
 278     TTree *treeInSignal = (TTree*) gDirectory->Get("TreeR");
 279     treeInSignal -> SetName("TreeR_InSignal");
 280     iPlane=0;
 281     while (treeInSignal->GetBranch(Form("Plane_%02d",iPlane))) {
 282       //        printf("Plane %02d\n",iPlane);
 283       treeInSignal ->SetBranchAddress(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneIn[iPlane]));
 284       TBranch *branch = treeOutTotal->GetBranch(Form("Plane_%02d",iPlane));
 285       if (!branch) treeOutTotal->Branch(Form("Plane_%02d",iPlane), &(fRecPointsPerPlaneOut[iPlane]));
 286       iPlane++;
 287     }
 288     iPlane=0;
 289     treeInSignal -> GetEntry(0);
 290     while (treeInSignal->GetBranch(Form("Plane_%02d",iPlane))) {
 291       Int_t nClusters = fRecPointsPerPlaneIn[iPlane]->GetEntries();
 292       for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
 293         //        printf("Cluster %4d\n",iCluster);
 294         AliMFTCluster *newCluster = (AliMFTCluster*) fRecPointsPerPlaneIn[iPlane]->At(iCluster);
 295         new ((*fRecPointsPerPlaneOut[iPlane])[fRecPointsPerPlaneOut[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
 296       }
 297       iPlane++;
 298     }
 299     for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 300       delete fRecPointsPerPlaneIn[jPlane];
 301     }
 302     delete treeInSignal;
 303     fileInSignal -> Close();
 304     delete fileInSignal;
 305
 306     if (eventExist) {
 307       treeOutTotal -> Fill();
 308       fileOutTotal -> mkdir(Form("Event%d",iEv));
 309       fileOutTotal -> cd(Form("Event%d",iEv));
 310       treeOutTotal -> SetName("TreeR");
 311       treeOutTotal -> Write();
 312       delete treeOutTotal;
 313       fileOutTotal -> Close();
 314       delete fileOutTotal;
 315       for (Int_t jPlane=0; jPlane<AliMFTConstants::fNMaxPlanes; jPlane++) {
 316         delete fRecPointsPerPlaneOut[jPlane];
 317       }
 318       iEv++;
 319     }
 320
 321   }
 322
 323 }
 324
 325 //====================================================================================================================================================
 326
 327 //  LocalWords:  nFiles