STEER/AliAODTagCreator.cxx

   1 /**************************************************************************
   2  * Author: Panos Christakoglou.                                           *
   3  * Contributors are mentioned in the code where appropriate.              *
   4  *                                                                        *
   5  * Permission to use, copy, modify and distribute this software and its   *
   6  * documentation strictly for non-commercial purposes is hereby granted   *
   7  * without fee, provided that the above copyright notice appears in all   *
   8  * copies and that both the copyright notice and this permission notice   *
   9  * appear in the supporting documentation. The authors make no claims     *
  10  * about the suitability of this software for any purpose. It is          *
  11  * provided "as is" without express or implied warranty.                  *
  12  **************************************************************************/
  13
  14 /* $Id$ */
  15
  16 //-----------------------------------------------------------------
  17 //           AliAODTagCreator class
  18 //   This is the class to deal with the tag creation (post process)
  19 //   Origin: Panos Christakoglou, UOA-CERN, Panos.Christakoglou@cern.ch
  20 //-----------------------------------------------------------------
  21
  22 //ROOT
  23 #include <Riostream.h>
  24 #include <TFile.h>
  25 #include <TString.h>
  26 #include <TTree.h>
  27 #include <TSystem.h>
  28 #include <TChain.h>
  29 #include <TLorentzVector.h>
  30 #include <TRefArray.h>
  31
  32 //ROOT-AliEn
  33 #include <TGrid.h>
  34 #include <TGridResult.h>
  35
  36 //AliRoot
  37 #include "AliRunTag.h"
  38 #include "AliEventTag.h"
  39 #include "AliFileTag.h"
  40 #include "AliPID.h"
  41 #include "AliAODEvent.h"
  42 #include "AliAODVertex.h"
  43 #include "AliLog.h"
  44
  45 #include "AliAODTagCreator.h"
  46
  47
  48 ClassImp(AliAODTagCreator)
  49
  50
  51 //______________________________________________________________________________
  52   AliAODTagCreator::AliAODTagCreator() :
  53     AliTagCreator(),
  54     fChain(0),
  55     fAODEvent(0),
  56     fTreeT(0),
  57     fRunTag(0),
  58     fTreeTEsd(0),
  59     fRunTagEsd(0)
  60 {
  61   //==============Default constructor for a AliAODTagCreator================
  62 }
  63
  64 //______________________________________________________________________________
  65 AliAODTagCreator::~AliAODTagCreator() {
  66 //================Default destructor for a AliAODTagCreator===================
  67   delete fChain;
  68   delete fAODEvent;
  69 }
  70
  71 //______________________________________________________________________________
  72 Bool_t AliAODTagCreator::ReadGridCollection(TGridResult *fresult) {
  73   // Reads the entry of the TGridResult and creates the tags
  74   Int_t nEntries = fresult->GetEntries();
  75
  76   TString alienUrl;
  77   const char* guid;
  78   const char* md5;
  79   const char* turl;
  80   Long64_t size = -1;
  81
  82   fChain = new TChain("aodTree");
  83
  84   for(Int_t i = 0; i < nEntries; i++) {
  85     alienUrl = fresult->GetKey(i,"turl");
  86     guid = fresult->GetKey(i,"guid");
  87     if(fresult->GetKey(i,"size")) size = atol (fresult->GetKey(i,"size"));
  88     md5 = fresult->GetKey(i,"md5");
  89     turl = fresult->GetKey(i,"turl");
  90     if(md5 && !strlen(guid)) md5 = 0;
  91     if(guid && !strlen(guid)) guid = 0;
  92
  93     fChain->Add(alienUrl);
  94   }//grid result loop
  95
  96   AliInfo(Form("AOD chain created......."));
  97   AliInfo(Form("Chain entries: %lld",fChain->GetEntries()));
  98
  99   CreateTag(fChain, "grid");
 100
 101   return kTRUE;
 102 }
 103
 104 //______________________________________________________________________________
 105 Bool_t AliAODTagCreator::ReadLocalCollection(const char *localpath, const char* pattern) {
 106   // Checks the different subdirs of the given local path and in the
 107   // case where it finds an AliAOD.root file it creates the tags
 108
 109   void *dira =  gSystem->OpenDirectory(localpath);
 110   Char_t fPath[512];
 111   const char * dirname  = 0x0;
 112   const char * filename = 0x0;
 113
 114   fChain = new TChain("aodTree");
 115
 116   while((dirname = gSystem->GetDirEntry(dira))) {
 117     snprintf(fPath,512,"%s/%s",localpath,dirname);
 118     void *dirb =  gSystem->OpenDirectory(fPath);
 119     while((filename = gSystem->GetDirEntry(dirb))) {
 120         TString bstr = dirname;
 121         if(bstr.Contains("..")) continue;
 122         if(strstr(filename,pattern)) {
 123             TString aodFileName;
 124             aodFileName = fPath;
 125             aodFileName += "/";
 126             aodFileName += pattern;
 127             fChain->Add(aodFileName);
 128         } //pattern check
 129     } //child directory's entry loop
 130   } //parent directory's entry loop
 131
 132   AliInfo(Form("AOD chain created......."));
 133   AliInfo(Form("Chain entries: %lld",fChain->GetEntries()));
 134
 135   CreateTag(fChain, "local");
 136
 137   return kTRUE;
 138 }
 139
 140 //______________________________________________________________________________
 141 Bool_t AliAODTagCreator::ReadCAFCollection(const char *filename) {
 142   // Temporary solution for CAF: Takes as an input the ascii file that
 143   // lists the AODs stored in the SE of the CAF and creates the tags.
 144
 145   // Open the input stream
 146   ifstream in;
 147   in.open(filename);
 148
 149   TString aodfile;
 150
 151   fChain = new TChain("aodTree");
 152
 153   // Read the input list of files and add them to the chain
 154   while(in.good()) {
 155     in >> aodfile;
 156     if (!aodfile.Contains("root")) continue; // protection
 157     fChain->Add(aodfile);
 158   }
 159
 160   AliInfo(Form("AOD chain created......."));
 161   AliInfo(Form("Chain entries: %lld",fChain->GetEntries()));
 162
 163   CreateTag(fChain, "proof");
 164
 165   return kTRUE;
 166 }
 167
 168 //__________________________________________________________________________
 169 void AliAODTagCreator::CreateAODTags(Int_t fFirstEvent, Int_t fLastEvent, TList */*grpList*/) {
 170   // Creates tag files for AODs
 171
 172     AliInfo(Form("Creating the AOD tags......."));
 173
 174     TFile *file = TFile::Open("AliAOD.root");
 175     if (!file || !file->IsOpen()) {
 176         AliError(Form("opening failed"));
 177         delete file;
 178         return ;
 179     }
 180
 181     fChain = new TChain("aodTree");
 182     fChain->Add("AliAOD.root");
 183
 184     fAODEvent = new AliAODEvent();
 185     fAODEvent->ReadFromTree(fChain);
 186
 187     Int_t lastEvent = 0;
 188     if(fLastEvent == -1) lastEvent = (Int_t)fChain->GetEntries();
 189     else lastEvent = fLastEvent;
 190
 191     char fileName[256];
 192     snprintf(fileName, 256, "Run%d.Event%d_%d.AOD.tag.root",
 193             fAODEvent->GetRunNumber(), fFirstEvent, lastEvent );
 194     AliInfo(Form("writing tags to file %s", fileName));
 195     AliDebug(1, Form("writing tags to file %s", fileName));
 196
 197     TFile* ftag = TFile::Open(fileName, "recreate");
 198
 199     fRunTag = new AliRunTag();
 200     fTreeT = new TTree("T","A Tree with event tags");
 201     TBranch * btag = fTreeT->Branch("AliTAG", &fRunTag);
 202     btag->SetCompressionLevel(9);
 203
 204     CreateTags();
 205
 206     ftag->cd();
 207     fTreeT->Fill();
 208     fRunTag->Clear();
 209     fTreeT->Write();
 210     ftag->Close();
 211     file->cd();
 212     file->Close();
 213 }
 214
 215 //_____________________________________________________________________________
 216 void AliAODTagCreator::CreateTag(TChain* chain, const char *type) {
 217
 218     // Private method that creates tag files
 219     //
 220
 221     //reading the esd tag file
 222     fTreeTEsd = new TChain("T");
 223     const char * tagPattern = "ESD.tag";
 224     // Open the working directory
 225     void * dirp = gSystem->OpenDirectory(gSystem->pwd());
 226     const char * name = 0x0;
 227     // Add all files matching *pattern* to the chain
 228     while((name = gSystem->GetDirEntry(dirp))) {
 229         if (strstr(name,tagPattern)) fTreeTEsd->Add(name);
 230     }//directory loop
 231     AliInfo(Form("Chained tag files: %lld",fTreeTEsd->GetEntries()));
 232
 233     fChain = chain;
 234
 235     TString fSession = type;
 236     TString fguid, fmd5, fturl;
 237     fAODEvent = new AliAODEvent();
 238     fAODEvent->ReadFromTree(fChain);
 239
 240     Int_t firstEvent = 0;
 241
 242     TString localFileName = "Run";
 243     localFileName += fAODEvent->GetRunNumber();
 244     localFileName += ".Event";
 245     localFileName += firstEvent;
 246     localFileName += "_";
 247     localFileName += chain->GetEntries(); //localFileName += "."; localFileName += Counter;
 248     localFileName += ".AOD.tag.root";
 249
 250     TString fileName;
 251
 252     if(fStorage == 0) {
 253         fileName = localFileName.Data();
 254         AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
 255     }
 256     else if(fStorage == 1) {
 257         TString alienLocation = "/alien";
 258         alienLocation += gGrid->Pwd();
 259         alienLocation += fgridpath.Data();
 260         alienLocation += "/";
 261         alienLocation +=  localFileName;
 262         alienLocation += "?se=";
 263         alienLocation += fSE.Data();
 264         fileName = alienLocation.Data();
 265         AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
 266     }
 267
 268     TFile* ftag = TFile::Open(fileName, "recreate");
 269
 270     fRunTag = new AliRunTag();
 271     fTreeT  = new TTree("T", "A Tree with event tags");
 272     TBranch * btag = fTreeT->Branch("AliTAG", &fRunTag);
 273     btag->SetCompressionLevel(9);
 274
 275     // Access information from esd tag
 276     fRunTagEsd = new AliRunTag();
 277     fTreeTEsd->SetBranchAddress("AliTAG",&fRunTagEsd);
 278     // Creating new information of aod
 279     AliInfo(Form("Creating the AOD tags......."));
 280     CreateTags(type);
 281     ftag->cd();
 282     fRunTag->Clear();
 283     fTreeT->Write();
 284     ftag->Close();
 285 }
 286
 287 void AliAODTagCreator::CreateTags(const char* type)
 288 {
 289     // Event loop for tag creation
 290     TString fturl;
 291     TString fguid;
 292     Int_t   oldRun = -1;
 293     fChain->GetEntry(0);
 294     TFile *f = fChain->GetFile();
 295     TString ftempGuid = f->GetUUID().AsString();
 296     // Loop over events
 297     Int_t nEvents = fChain->GetEntries();
 298     Int_t ntags    = 0;
 299     Int_t tagentry = 0;
 300     //    const TClonesArray *evTagList = 0;
 301     TString foldguid = "";
 302
 303     for (Int_t iEventNumber = 0; iEventNumber < nEvents; iEventNumber++) {
 304         // Copy old tag information
 305         if (iEventNumber >= ntags) {
 306             fTreeTEsd->GetEntry(tagentry++);
 307             fRunTag->CopyStandardContent(fRunTagEsd);
 308 //          evTagList = fRunTagEsd->GetEventTags();
 309 //          ntags += evTagList->GetEntries();
 310             ntags = fRunTagEsd->GetNEvents();
 311         }
 312
 313         // Create a new Tag
 314         AliEventTag* evTag = new AliEventTag();
 315         // Read event
 316         fChain->GetEntry(iEventNumber);
 317         if (iEventNumber == 0) oldRun = fAODEvent->GetRunNumber();
 318         // Reference to the input file
 319         TFile *file = fChain->GetFile();
 320         //      const TUrl *url = file->GetEndpointUrl();
 321         fguid = file->GetUUID().AsString();
 322
 323 //      if (!strcmp(type,"grid")) {
 324 //          TString fturltemp = "alien://"; fturltemp += url->GetFile();
 325 //          fturl = fturltemp(0,fturltemp.Index(".root",5,0,TString::kExact)+5);
 326 //      } else {
 327 //          fturl = url->GetFile();
 328 //      }
 329         fturl = file->GetName();
 330
 331         fAODEvent->GetStdContent();
 332
 333         // Fill the event tag from the aod informatiom
 334         FillEventTag(fAODEvent, evTag);
 335         // Set the event and input file references
 336         //evTag->SetEventId(iEventNumber+1);
 337
 338         // **** FIXME ****
 339 //      evTag->SetGUID(fguid);
 340 //      if(!strcmp(type,"grid")) {
 341 //          evTag->SetMD5("");
 342 //          evTag->SetTURL(fturl);
 343 //          evTag->SetSize(0);
 344 //          }
 345 //      else evTag->SetPath(fturl);
 346         //  **** FIXME ****
 347
 348         // Check if a new run has to be created
 349         // File has changed
 350         if(fguid != ftempGuid) {
 351             ftempGuid = fguid;
 352             fTreeT->Fill();
 353             fRunTag->Clear("");
 354
 355             AliFileTag *nftag = new AliFileTag();
 356
 357             // if(fSession == "grid") {
 358               nftag->SetMD5("");
 359               nftag->SetTURL(fturl);
 360               nftag->SetSize(0);
 361             // }
 362             // else {
 363             //   nftag->SetPath(fturl);
 364             //   nftag->SetSize(0);
 365             //   nftag->SetMD5("");
 366             //   nftag->SetTURL(fturl);
 367             // }
 368
 369             if (fRunTag->GetFileId(fguid) > -1)
 370               AliFatal("Adding a file which is already in the RunTag.");
 371
 372             fRunTag->AddFileTag(nftag);
 373
 374         }
 375
 376         // Run# has changed
 377         if (oldRun != (fAODEvent->GetRunNumber()))
 378         {
 379             oldRun = fAODEvent->GetRunNumber();
 380
 381             fTreeT->Fill();
 382             fRunTag->Clear("");
 383             ftempGuid = fguid;
 384             fTreeT->Fill();
 385             fRunTag->Clear("");
 386
 387             AliFileTag *nftag = new AliFileTag();
 388
 389             // if(fSession == "grid") {
 390               nftag->SetMD5("");
 391               nftag->SetTURL(fturl);
 392               nftag->SetSize(0);
 393             // }
 394             // else {
 395             //   nftag->SetPath(fturl);
 396             //   nftag->SetSize(0);
 397             //   nftag->SetMD5("");
 398             //   nftag->SetTURL(fturl);
 399             // }
 400
 401             if (fRunTag->GetFileId(fguid) > -1)
 402               AliFatal("Adding a file which is already in the RunTag.");
 403
 404             fRunTag->AddFileTag(nftag);
 405
 406         }
 407
 408         // Add the event tag
 409         fRunTag->AddEventTag(*evTag);
 410         delete evTag;
 411         // Last event
 412         if(iEventNumber+1 == fChain->GetEntries()) {
 413             fTreeT->Fill();
 414             fRunTag->Clear("");
 415         }
 416     }//event loop
 417 }
 418
 419
 420
 421 void AliAODTagCreator::FillEventTag(AliAODEvent* aod, AliEventTag* evTag)
 422 {
 423 //
 424 // Fill the event tag information
 425     //
 426     fAODEvent = aod;
 427
 428     //
 429     Float_t fLowPtCut      =  1.0;
 430     Float_t fHighPtCut     =  3.0;
 431     Float_t fVeryHighPtCut = 10.0;
 432     ////////////
 433     Double_t partFrac[10] = {0.01, 0.01, 0.85, 0.10, 0.05, 0., 0., 0., 0., 0.};
 434
 435     // Creates the tags for all the events in a given AOD file
 436     Int_t   ntrack = 0;
 437     Int_t   nPos = 0, nNeg = 0, nNeutr =0;
 438     Int_t   nKinks = 0, nV0s = 0, nCascades = 0;
 439     Int_t   nK0s = 0, nNeutrons = 0, nPi0s = 0, nGamas = 0;
 440     Int_t   nProtons = 0,  nKaons = 0, nPions = 0, nMuons = 0, nElectrons = 0, nFWMuons = 0;
 441     Int_t   nCh1GeV = 0, nCh3GeV = 0, nCh10GeV = 0;
 442     Int_t   nMu1GeV = 0, nMu3GeV = 0, nMu10GeV = 0;
 443     Int_t   nEl1GeV = 0, nEl3GeV = 0, nEl10GeV = 0;
 444     Float_t maxPt =  .0, etamaxPt = -999., phimaxPt = -999., meanPt = .0, totalP =  .0;
 445
 446     TRefArray tmp;
 447
 448
 449     // Primary Vertex
 450     AliAODVertex *pVertex = fAODEvent->GetPrimaryVertex();
 451     if (pVertex) {
 452         evTag->SetVertexX(pVertex->GetX());
 453         evTag->SetVertexY(pVertex->GetY());
 454         evTag->SetVertexZ(pVertex->GetZ());
 455         Double_t covmatrix[6];
 456         pVertex->GetCovarianceMatrix(covmatrix);
 457         evTag->SetVertexZError(sqrt(covmatrix[5]));
 458     }
 459     // loop over vertices
 460     Int_t nVtxs = fAODEvent->GetNumberOfVertices();
 461     for (Int_t nVtx = 0; nVtx < nVtxs; nVtx++) {
 462         AliAODVertex *vertex = fAODEvent->GetVertex(nVtx);
 463         if(vertex->GetType() == 1) nKinks    += 1;
 464         if(vertex->GetType() == 2) nV0s      += 1;
 465         if(vertex->GetType() == 3) nCascades += 1;
 466     }
 467     Int_t nTracks = fAODEvent->GetNTracks();
 468     for (Int_t nTr = 0; nTr < nTracks; nTr++) {
 469         AliAODTrack *track = fAODEvent->GetTrack(nTr);
 470
 471         Double_t fPt = track->Pt();
 472         if(fPt > maxPt) {
 473             maxPt = fPt;
 474             etamaxPt = track->Eta();
 475             phimaxPt = track->Phi();
 476         }
 477
 478         if(track->Charge() > 0) {
 479             nPos++;
 480             if(fPt > fLowPtCut)      nCh1GeV++;
 481             if(fPt > fHighPtCut)     nCh3GeV++;
 482             if(fPt > fVeryHighPtCut) nCh10GeV++;
 483         }
 484         if(track->Charge() < 0) {
 485             nNeg++;
 486             if(fPt > fLowPtCut)      nCh1GeV++;
 487             if(fPt > fHighPtCut)     nCh3GeV++;
 488             if(fPt > fVeryHighPtCut) nCh10GeV++;
 489             }
 490         if(track->Charge() == 0) nNeutr++;
 491         //PID
 492         const Double32_t *prob = track->PID();
 493         Double_t rcc = 0.0;
 494         for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
 495         if(rcc == 0.0) continue;
 496         //Bayes' formula
 497         Double_t w[10];
 498         for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
 499
 500         //protons
 501         if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
 502         //kaons
 503         if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
 504         //pions
 505         if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
 506         //muons
 507         if ((w[1]>w[4])&&(w[1]>w[3])&&(w[1]>w[2])&&(w[1]>w[0])) {
 508             nMuons++;
 509             if(fPt > fLowPtCut)      nMu1GeV++;
 510             if(fPt > fHighPtCut)     nMu3GeV++;
 511             if(fPt > fVeryHighPtCut) nMu10GeV++;
 512         }
 513         //electrons
 514         if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
 515             nElectrons++;
 516             if(fPt > fLowPtCut) nEl1GeV++;
 517             if(fPt > fHighPtCut) nEl3GeV++;
 518             if(fPt > fVeryHighPtCut) nEl10GeV++;
 519         }
 520         totalP += track->P();
 521         meanPt += fPt;
 522         ntrack++;
 523         // forward muons (in the dimuon spectrometer)
 524         if(track->IsMuonTrack()) nFWMuons++;
 525
 526     }//track loop
 527     //
 528     // Fill the event tags
 529     if(ntrack != 0)
 530         meanPt = meanPt/ntrack;
 531
 532
 533     evTag->SetNumOfTracks(nTracks);
 534     evTag->SetNumOfPosTracks(nPos);
 535     evTag->SetNumOfNegTracks(nNeg);
 536     evTag->SetNumOfNeutrTracks(nNeutr);
 537
 538     evTag->SetNumOfV0s(nV0s);
 539     evTag->SetNumOfCascades(nCascades);
 540     evTag->SetNumOfKinks(nKinks);
 541
 542     evTag->SetNumOfProtons(nProtons);
 543     evTag->SetNumOfKaons(nKaons);
 544     evTag->SetNumOfPions(nPions);
 545     evTag->SetNumOfMuons(nMuons);
 546     evTag->SetNumOfFWMuons(nFWMuons);
 547     evTag->SetNumOfElectrons(nElectrons);
 548     evTag->SetNumOfPhotons(nGamas);
 549     evTag->SetNumOfPi0s(nPi0s);
 550     evTag->SetNumOfNeutrons(nNeutrons);
 551     evTag->SetNumOfKaon0s(nK0s);
 552
 553     evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
 554     evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
 555     evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
 556     evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
 557     evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
 558     evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
 559     evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
 560     evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
 561     evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
 562
 563     tmp.Clear();
 564     evTag->SetNumOfPHOSClusters(fAODEvent->GetPHOSClusters(&tmp));
 565     tmp.Clear();
 566     evTag->SetNumOfEMCALClusters(fAODEvent->GetEMCALClusters(&tmp));
 567
 568     evTag->SetTotalMomentum(totalP);
 569     evTag->SetMeanPt(meanPt);
 570     evTag->SetMaxPt(maxPt);
 571     evTag->SetEtaMaxPt(etamaxPt);
 572     evTag->SetPhiMaxPt(phimaxPt);
 573 }