TRD/Macros/AliTRDCreateVector2DSum.C

   1 #if !defined( __CINT__) || defined(__MAKECINT__)
   2
   3 #include <vector>
   4 #include <TChain.h>
   5 #include <TObjArray.h>
   6 #include <TFile.h>
   7 #include <TGraphErrors.h>
   8 #include <TTree.h>
   9 #include <TDirectory.h>
  10 #include <TROOT.h>
  11 #include <TString.h>
  12 #include <Riostream.h>
  13 #include <TSystem.h>
  14 #include <TH1F.h>
  15 #include <AliCDBManager.h>
  16 #include "../TRD/AliTRDCalibra.h"
  17
  18 #endif
  19
  20
  21 Bool_t AliTRDCreateVector2DSum(const char* variablecali, const char* noime, const char* dire, const char* namefile){
  22   //
  23   // After having simulated and reconstructed events in subrepertories 000%d of dire
  24   // this macro searchs in the subdirectories the file TRD.calibration.root
  25   // takes the vectors and merge them
  26   // variablecali can be treeCH2d, treePH2d or treePRF2d
  27   // noime is the name of the "sum" tree
  28   // namefile is the name of the file where noime will be written
  29   //
  30
  31
  32   // The tree
  33   TChain *treeChain = new TChain(variablecali);
  34
  35
  36   //Variables
  37   TObjArray *vectorplace = new TObjArray();
  38   TObjArray *whereinthechain = new TObjArray();
  39
  40
  41   //TRDCalibra
  42   AliCDBManager *man = AliCDBManager::Instance();
  43   man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
  44   man->SetRun(0);
  45   AliTRDCalibra *calibra = AliTRDCalibra::Instance();
  46
  47   //The patterns
  48   const char *name="TRD.calibration.root";
  49   const char *patterndir="0";
  50   const char *namesubdir = 0x0;
  51   char dir[200];
  52
  53   //Open the directory dire
  54   void * dircu = gSystem->OpenDirectory(dire);
  55   while((namesubdir = gSystem->GetDirEntry(dircu))) {
  56     if(strstr(namesubdir,patterndir)) {
  57       sprintf(dir,"%s/%s",dire,namesubdir);
  58       printf("process subdirectories: %s\n",dir);
  59
  60       //full name of the file and tree
  61       char fullname[255] = "";
  62
  63       sprintf(fullname,"%s/%s",dir,name);
  64       printf("Process file: %s\n",fullname);
  65       TFile *file = new TFile(fullname,"READ");
  66       if(!file) continue;
  67       TTree *treecurrent = (TTree *) file->Get(variablecali);
  68       if(!treecurrent) continue;
  69       gDirectory = gROOT;
  70       TObjArray *vectorcurrent = calibra->ConvertTreeVector(treecurrent);
  71       printf("Size of the current tree: %d\n",(Int_t) vectorcurrent->GetEntriesFast());
  72       printf("Size of whereinthechain: %d\n",(Int_t) whereinthechain->GetEntriesFast());
  73       printf("Size of the chain: %d\n", (Int_t) treeChain->GetEntries());
  74       Int_t j = (Int_t) treeChain->GetEntries();
  75         for(Int_t jui = 0; jui < (Int_t) vectorcurrent->GetEntriesFast(); jui++){
  76           //Search if already found
  77           Int_t place = calibra->SearchInTreeVector(vectorplace,((AliTRDCalibra::AliTRDPlace *)vectorcurrent->At(jui))->GetPlace());
  78           //Create a new element in the two std vectors
  79           if(place == -1){
  80             AliTRDCalibra::AliTRDPlace *placejui = new AliTRDCalibra::AliTRDPlace();
  81             placejui->SetPlace(j+jui);
  82             TObjArray *chainplace = new TObjArray();
  83             chainplace->Add((TObject *) placejui);
  84             vectorplace->Add((TObject *) (vectorcurrent->At(jui)));
  85             whereinthechain->Add((TObject *) chainplace);
  86           }
  87           //Update the element at the place "place" in the std vector whereinthechain
  88           else {
  89             AliTRDCalibra::AliTRDPlace *placejui = new AliTRDCalibra::AliTRDPlace();
  90             placejui->SetPlace((j+jui));
  91             TObjArray *chainplace = ((TObjArray *) whereinthechain->At(place));
  92             chainplace->Add((TObject *) placejui);
  93             whereinthechain->AddAt((TObject *)chainplace, place);
  94           }
  95         }
  96         treeChain->AddFile(fullname);
  97         delete file;
  98     }//if pattern
  99   }// in the sub
 100
 101   //Take care of the profile
 102   const char* pattern = "P";
 103   TTree *tree;
 104
 105   if(!strstr(variablecali,pattern)){
 106     //Ready to read the chain
 107     TH1F *his = 0x0;
 108     treeChain->SetBranchAddress("histo",&his);
 109
 110     //Initialise the final tree
 111     Int_t group = -1;
 112     TH1F *histsum = 0x0;
 113
 114     tree = new TTree(noime,noime);
 115     tree->Branch("groupnumber",&group,"groupnumber/I");
 116     tree->Branch("histo","TH1F",&histsum,32000,0);
 117
 118     //Init histsum
 119     if(treeChain->GetEntries() < 1) return kFALSE;
 120
 121     printf("FINAL Size of the chain: %d\n", (Int_t) treeChain->GetEntriesFast());
 122     printf("FINAL Size of vectorplace: %d\n", (Int_t) vectorplace->GetEntriesFast());
 123     for(Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++){
 124       group = ((AliTRDCalibra::AliTRDPlace *)vectorplace->At(h))->GetPlace();
 125       TObjArray *chainplace = ((TObjArray *)whereinthechain->At(h));
 126       treeChain->GetEntry(((AliTRDCalibra::AliTRDPlace *)chainplace->At(0))->GetPlace());
 127       //Init for the first time
 128       if(h == 0)  {
 129         histsum = new TH1F("","",his->GetXaxis()->GetNbins(),his->GetXaxis()->GetBinLowEdge(1),his->GetXaxis()->GetBinUpEdge(his->GetXaxis()->GetNbins()));
 130         histsum->Sumw2();
 131       }
 132       //Reset for each new group
 133       histsum->SetEntries(0);
 134       for(Int_t l = 0; l <= histsum->GetXaxis()->GetNbins(); l++){
 135         histsum->SetBinContent(l,0.0);
 136         histsum->SetBinError(l,0.0);
 137       }
 138       histsum->Add(his,1);
 139       if((Int_t) chainplace->GetEntriesFast() > 1){
 140         for(Int_t s = 1; s < (Int_t)chainplace->GetEntriesFast(); s++){
 141           treeChain->GetEntry(((AliTRDCalibra::AliTRDPlace *)chainplace->At(s))->GetPlace());
 142           histsum->Add(his,1);
 143         }
 144       }
 145       tree->Fill();
 146     }
 147   }
 148   else {
 149     //Ready to read the chain
 150     TGraphErrors *his = 0x0;
 151     treeChain->SetBranchAddress("histo",&his);
 152
 153     //Initialise the final tree
 154     Int_t group = -1;
 155     TGraphErrors *histsum = 0x0;
 156     Double_t *xref = 0x0;
 157
 158
 159     tree = new TTree(noime,noime);
 160     tree->Branch("groupnumber",&group,"groupnumber/I");
 161     tree->Branch("histo","TGraphErrors",&histsum,32000,0);
 162
 163     //Init histsum
 164     if(treeChain->GetEntries() < 1) return kFALSE;
 165
 166     printf("FINAL Size of the chain: %d\n", (Int_t) treeChain->GetEntriesFast());
 167     printf("FINAL Size of vectorplace: %d\n", (Int_t) vectorplace->GetEntriesFast());
 168     for(Int_t h = 0; h < (Int_t) vectorplace->GetEntriesFast(); h++){
 169       group = ((AliTRDCalibra::AliTRDPlace *)vectorplace->At(h))->GetPlace();
 170       TObjArray *chainplace = ((TObjArray *)whereinthechain->At(h));
 171       treeChain->GetEntry(((AliTRDCalibra::AliTRDPlace *)chainplace->At(0))->GetPlace());
 172       //Init for the fisrt time
 173       Int_t nbins = his->GetN();
 174       Double_t *x;
 175       x = new Double_t[nbins];
 176       xref = his->GetX();
 177       Double_t *ex;
 178       ex = new Double_t[nbins];
 179       Double_t *y;
 180       y = new Double_t[nbins];
 181       Double_t *ey;
 182       ey = new Double_t[nbins];
 183
 184       for(Int_t lo = 0; lo < nbins; lo++){
 185         x[lo] = xref[lo];
 186         ex[lo] = 0.0;
 187         y[lo]  = 0.0;
 188         ey[lo] = 0.0;
 189       }
 190       delete histsum;
 191       histsum = new TGraphErrors(nbins,x,y,ex,ey);
 192       //Reset for each group
 193       Double_t *xp;
 194       xp = histsum->GetX();
 195       Double_t *exp;
 196       exp = histsum->GetEX();
 197       Double_t *yp;
 198       yp = histsum->GetY();
 199       Double_t *eyp;
 200       eyp = histsum->GetEY();
 201
 202       //Add the first
 203       histsum = calibra->AddProfiles(his,histsum);
 204       if((Int_t) chainplace->GetEntriesFast() > 1){
 205         for(Int_t s = 1; s < (Int_t)chainplace->GetEntriesFast(); s++){
 206           treeChain->GetEntry(((AliTRDCalibra::AliTRDPlace *)chainplace->At(s))->GetPlace());
 207           histsum = calibra->AddProfiles(his,histsum);
 208         }
 209       }
 210       tree->Fill();
 211     }
 212
 213   }
 214
 215
 216   //Write in a file
 217   TFile *fout = TFile::Open((const char*) namefile,"UPDATE");
 218   fout->WriteTObject(tree,tree->GetName(),(Option_t *) "kOverWrite");
 219   fout->Close();
 220
 221   return kTRUE;
 222 }