[u/mrichter/AliRoot.git] / STEER / AliTagCreator.cxx

/**************************************************************************
 * Author: Panos Christakoglou.                                           *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

//-----------------------------------------------------------------
//           AliTagCreator class
//   This is the class to deal with the tag creation (post process)
//   Origin: Panos Christakoglou, UOA-CERN, Panos.Christakoglou@cern.ch
//-----------------------------------------------------------------

//ROOT
#include <TFile.h>
#include <TString.h>
#include <TTree.h>
#include <TSystem.h>
#include <TChain.h>
#include <TLorentzVector.h>

//ROOT-AliEn
#include <TGrid.h>
#include <TGridResult.h>

//AliRoot
#include "AliRunTag.h"
#include "AliEventTag.h"
#include "AliESD.h"
#include "AliESDVertex.h"
#include "AliLog.h"


#include "AliTagCreator.h"


ClassImp(AliTagCreator)


//______________________________________________________________________________
AliTagCreator::AliTagCreator() //local mode
{
  //==============Default constructor for a AliTagCreator==================
  fgridpath = "";
  fSE = "";   
  fStorage = 0; 
}

//______________________________________________________________________________
AliTagCreator::~AliTagCreator()
{
//================Default destructor for a AliTagCreator=======================
}

//______________________________________________________________________________
void AliTagCreator::SetStorage(Int_t storage)
{
  // Sets correctly the storage: 0 for local, 1 for GRID
  fStorage = storage;
  if(fStorage == 0)
    AliInfo(Form("Tags will be stored locally...."));
  if(fStorage == 1)
    AliInfo(Form("Tags will be stored in the grid...."));
  if((fStorage != 0)&&(fStorage != 1))
    {
      AliInfo(Form("Storage was not properly set!!!"));
      abort();
    }  
}


//______________________________________________________________________________
Bool_t AliTagCreator::ReadESDCollection(TGridResult *fresult)
{
  // Reads the entry of the TGridResult and creates the tags
  Int_t nEntries = fresult->GetEntries();

  TString alienUrl;
  const char *guid;
  const char *md5;
  const char *turl;
  Long64_t size = -1;

  Int_t counter = 0;
  for(Int_t i = 0; i < nEntries; i++)
    {
      alienUrl = fresult->GetKey(i,"turl");
      guid = fresult->GetKey(i,"guid");
      if(fresult->GetKey(i,"size"))
	size = atol (fresult->GetKey(i,"size"));
      md5 = fresult->GetKey(i,"md5");
      turl = fresult->GetKey(i,"turl");
      if(md5 && !strlen(guid))
	md5 = 0;
      if(guid && !strlen(guid))
	guid = 0;

      TFile *f = TFile::Open(alienUrl,"READ");
      CreateTag(f,guid,md5,turl,size,counter);
      f->Close();
      delete f;	 
      counter += 1;
    }//grid result loop

  return kTRUE;
}

//__________________________________________________________________________
Bool_t AliTagCreator::MergeTags()
{
  //Merges the tags and stores the merged tag file 
  //locally if fStorage=0 or in the grid if fStorage=1
  AliInfo(Form("Merging tags....."));
  TChain *fgChain = new TChain("T");

  if(fStorage == 0) {
    const char * tagPattern = "tag";
    // Open the working directory
    void * dirp = gSystem->OpenDirectory(gSystem->pwd());
    const char * name = 0x0;
    // Add all files matching *pattern* to the chain
    while((name = gSystem->GetDirEntry(dirp))) {
      if (strstr(name,tagPattern))       
	fgChain->Add(name);  
    }//directory loop
    AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
  }//local mode

  else if(fStorage == 1) {
    TString alienLocation = gGrid->Pwd();
    alienLocation += fgridpath.Data();
    alienLocation += "/";

    TGridResult *tagresult = gGrid->Query(alienLocation,"*tag.root","","");
    Int_t nEntries = tagresult->GetEntries();
    for(Int_t i = 0; i < nEntries; i++) {
      TString alienUrl = tagresult->GetKey(i,"turl");
      fgChain->Add(alienUrl);
    }//grid result loop      
    AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
  }//grid mode
 
  AliRunTag *tag = new AliRunTag;
  AliEventTag *evTag = new AliEventTag;
  fgChain->SetBranchAddress("AliTAG",&tag);
   
  //Defining new tag objects
  AliRunTag *newTag = new AliRunTag();
  TTree ttag("T","A Tree with event tags");
  TBranch * btag = ttag.Branch("AliTAG", &newTag);
  btag->SetCompressionLevel(9);
  for(Int_t iTagFiles = 0; iTagFiles < fgChain->GetEntries(); iTagFiles++) {
    fgChain->GetEntry(iTagFiles);
    newTag->SetRunId(tag->GetRunId());
    const TClonesArray *tagList = tag->GetEventTags();
    for(Int_t j = 0; j < tagList->GetEntries(); j++) {
      evTag = (AliEventTag *) tagList->At(j);
      newTag->AddEventTag(*evTag);
    }
    ttag.Fill();
    newTag->Clear();
  }//tag file loop 
  
  TString localFileName = "Run"; localFileName += tag->GetRunId(); 
  localFileName += ".Merged"; localFileName += ".ESD.tag.root";
     
  TString alienFileName = "/alien";
  alienFileName += gGrid->Pwd();
  alienFileName += fgridpath.Data();
  alienFileName += "/";
  alienFileName +=  localFileName;
  alienFileName += "?se=";
  alienFileName += fSE.Data();

  TString filename = 0x0;
  
  if(fStorage == 0) {
    filename = localFileName.Data();      
    AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
  } 
  if(fStorage == 1) {
    filename = alienFileName.Data();
    AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));     
  }
  
  TFile* ftag = TFile::Open(filename, "recreate");
  ftag->cd();
  ttag.Write();
  ftag->Close();

  delete tag;
  delete newTag;

  return kTRUE;
}

//_____________________________________________________________________________
void AliTagCreator::CreateTag(TFile* file, const char *guid, const char *md5, const char *turl, Long64_t size, Int_t Counter)
{
  /////////////
  //muon code//
  ////////////
  Double_t fMUONMASS = 0.105658369;
  //Variables
  Double_t fX,fY,fZ ;
  Double_t fThetaX, fThetaY, fPyz, fChisquare;
  Double_t fPxRec, fPyRec, fPzRec, fEnergy;
  Int_t fCharge;
  TLorentzVector fEPvector;

  Float_t fZVertexCut = 10.0; 
  Float_t fRhoVertexCut = 2.0; 

  Float_t fLowPtCut = 1.0;
  Float_t fHighPtCut = 3.0;
  Float_t fVeryHighPtCut = 10.0;
  ////////////

  Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};

  // Creates the tags for all the events in a given ESD file
  Int_t ntrack;
  Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
  Int_t nPos, nNeg, nNeutr;
  Int_t nK0s, nNeutrons, nPi0s, nGamas;
  Int_t nCh1GeV, nCh3GeV, nCh10GeV;
  Int_t nMu1GeV, nMu3GeV, nMu10GeV;
  Int_t nEl1GeV, nEl3GeV, nEl10GeV;
  Float_t maxPt = .0, meanPt = .0, totalP = .0;
  Int_t fVertexflag;
  TString fVertexName;

  AliRunTag *tag = new AliRunTag();
  AliEventTag *evTag = new AliEventTag();
  TTree ttag("T","A Tree with event tags");
  TBranch * btag = ttag.Branch("AliTAG", &tag);
  btag->SetCompressionLevel(9);
  
  AliInfo(Form("Creating the tags......."));	
  
  Int_t firstEvent = 0,lastEvent = 0;
  TTree *t = (TTree*) file->Get("esdTree");
  TBranch * b = t->GetBranch("ESD");
  AliESD *esd = 0;
  b->SetAddress(&esd);
  
  tag->SetRunId(esd->GetRunNumber());
  
  Int_t iNumberOfEvents = b->GetEntries();
  for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
    ntrack = 0;
    nPos = 0;
    nNeg = 0;
    nNeutr =0;
    nK0s = 0;
    nNeutrons = 0;
    nPi0s = 0;
    nGamas = 0;
    nProtons = 0;
    nKaons = 0;
    nPions = 0;
    nMuons = 0;
    nElectrons = 0;	  
    nCh1GeV = 0;
    nCh3GeV = 0;
    nCh10GeV = 0;
    nMu1GeV = 0;
    nMu3GeV = 0;
    nMu10GeV = 0;
    nEl1GeV = 0;
    nEl3GeV = 0;
    nEl10GeV = 0;
    maxPt = .0;
    meanPt = .0;
    totalP = .0;
    fVertexflag = 1;
    
    b->GetEntry(iEventNumber);
    const AliESDVertex * vertexIn = esd->GetVertex();
    fVertexName = vertexIn->GetName();
    if(fVertexName == "default") fVertexflag = 0;

    for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
      AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
      UInt_t status = esdTrack->GetStatus();
      
      //select only tracks with ITS refit
      if ((status&AliESDtrack::kITSrefit)==0) continue;
      //select only tracks with TPC refit
      if ((status&AliESDtrack::kTPCrefit)==0) continue;
      
      //select only tracks with the "combined PID"
      if ((status&AliESDtrack::kESDpid)==0) continue;
      Double_t p[3];
      esdTrack->GetPxPyPz(p);
      Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
      Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
      totalP += momentum;
      meanPt += fPt;
      if(fPt > maxPt) maxPt = fPt;
      
      if(esdTrack->GetSign() > 0) {
	nPos++;
	if(fPt > fLowPtCut) nCh1GeV++;
	if(fPt > fHighPtCut) nCh3GeV++;
	if(fPt > fVeryHighPtCut) nCh10GeV++;
      }
      if(esdTrack->GetSign() < 0) {
	nNeg++;
	if(fPt > fLowPtCut) nCh1GeV++;
	if(fPt > fHighPtCut) nCh3GeV++;
	if(fPt > fVeryHighPtCut) nCh10GeV++;
      }
      if(esdTrack->GetSign() == 0) nNeutr++;
      
      //PID
      Double_t prob[5];
      esdTrack->GetESDpid(prob);
      
      Double_t rcc = 0.0;
      for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
      if(rcc == 0.0) continue;
      //Bayes' formula
      Double_t w[5];
      for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
      
      //protons
      if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
      //kaons
      if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
      //pions
      if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++; 
      //electrons
      if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
	nElectrons++;
	if(fPt > fLowPtCut) nEl1GeV++;
	if(fPt > fHighPtCut) nEl3GeV++;
	if(fPt > fVeryHighPtCut) nEl10GeV++;
      }	  
      ntrack++;
    }//esd track loop
    
    /////////////
    //muon code//
    ////////////
    Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
    // loop over all reconstructed tracks (also first track of combination)
    for (Int_t iTrack = 0; iTrack <  nMuonTracks;  iTrack++) {
      AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
      if (muonTrack == 0x0) continue;
      
      // Coordinates at vertex
      fZ = muonTrack->GetZ(); 
      fY = muonTrack->GetBendingCoor();
      fX = muonTrack->GetNonBendingCoor(); 
      
      fThetaX = muonTrack->GetThetaX();
      fThetaY = muonTrack->GetThetaY();
      
      fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
      fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
      fPxRec = fPzRec * TMath::Tan(fThetaX);
      fPyRec = fPzRec * TMath::Tan(fThetaY);
      fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
      
      //ChiSquare of the track if needed
      fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
      fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
      fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
      
      // total number of muons inside a vertex cut 
      if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
	nMuons++;
	if(fEPvector.Pt() > fLowPtCut) {
	  nMu1GeV++; 
	  if(fEPvector.Pt() > fHighPtCut) {
	    nMu3GeV++; 
	    if (fEPvector.Pt() > fVeryHighPtCut) {
	      nMu10GeV++;
	    }
	  }
	}
      }
    }//muon track loop
    
    // Fill the event tags 
    if(ntrack != 0)
      meanPt = meanPt/ntrack;
    
    evTag->SetEventId(iEventNumber+1);
    evTag->SetGUID(guid);
    evTag->SetMD5(md5);
    evTag->SetTURL(turl);
    evTag->SetSize(size);
    evTag->SetVertexX(vertexIn->GetXv());
    evTag->SetVertexY(vertexIn->GetYv());
    evTag->SetVertexZ(vertexIn->GetZv());
    evTag->SetVertexZError(vertexIn->GetZRes());
    evTag->SetVertexFlag(fVertexflag);
    
    evTag->SetT0VertexZ(esd->GetT0zVertex());
    
    evTag->SetTrigger(esd->GetTriggerMask());
    
    evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
    evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
    evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
    evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
    evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
    evTag->SetNumOfParticipants(esd->GetZDCParticipants());
    
    
    evTag->SetNumOfTracks(esd->GetNumberOfTracks());
    evTag->SetNumOfPosTracks(nPos);
    evTag->SetNumOfNegTracks(nNeg);
    evTag->SetNumOfNeutrTracks(nNeutr);
    
    evTag->SetNumOfV0s(esd->GetNumberOfV0s());
    evTag->SetNumOfCascades(esd->GetNumberOfCascades());
    evTag->SetNumOfKinks(esd->GetNumberOfKinks());
    evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
    
    evTag->SetNumOfProtons(nProtons);
    evTag->SetNumOfKaons(nKaons);
    evTag->SetNumOfPions(nPions);
    evTag->SetNumOfMuons(nMuons);
    evTag->SetNumOfElectrons(nElectrons);
    evTag->SetNumOfPhotons(nGamas);
    evTag->SetNumOfPi0s(nPi0s);
    evTag->SetNumOfNeutrons(nNeutrons);
    evTag->SetNumOfKaon0s(nK0s);
    
    evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
    evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
    evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
    evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
    evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
    evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
    evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
    evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
    evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
    
    evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
    evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
    
    evTag->SetTotalMomentum(totalP);
    evTag->SetMeanPt(meanPt);
    evTag->SetMaxPt(maxPt);
    
    tag->AddEventTag(*evTag);
  }//event loop
  lastEvent = iNumberOfEvents;
  
  t->Delete("");
  
  ttag.Fill();
  tag->Clear();

  TString localFileName = "Run"; localFileName += tag->GetRunId(); 
  localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
  localFileName += ".ESD.tag.root";

  TString alienLocation = "/alien";
  alienLocation += gGrid->Pwd();
  alienLocation += fgridpath.Data();
  alienLocation += "/";
  alienLocation +=  localFileName;
  alienLocation += "?se=";
  alienLocation += fSE.Data();

  TString fileName;
  
  if(fStorage == 0)
    {
      fileName = localFileName.Data();      
      AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
    }
  if(fStorage == 1)
    {
      fileName = alienLocation.Data();
      AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
    }

  TFile* ftag = TFile::Open(fileName, "recreate");
  ftag->cd();
  ttag.Write();
  ftag->Close();

  delete ftag;
  delete esd;

  delete tag;
}
Commit	Line	Data
b45e5084	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	// AliTagCreator 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
b45e5084	23	#include <TFile.h>
b45e5084	24	#include <TString.h>
cb1645b7	25	#include <TTree.h>
4fd84456	26	#include <TSystem.h>
4fd84456	27	#include <TChain.h>
2bdb9d38	28	#include <TLorentzVector.h>
b45e5084	29
	30	//ROOT-AliEn
	31	#include <TGrid.h>
b45e5084	32	#include <TGridResult.h>
b45e5084	33
	34	//AliRoot
	35	#include "AliRunTag.h"
	36	#include "AliEventTag.h"
	37	#include "AliESD.h"
	38	#include "AliESDVertex.h"
b45e5084	39	#include "AliLog.h"
	40
	41
	42	#include "AliTagCreator.h"
	43
	44
	45	ClassImp(AliTagCreator)
	46
	47
	48	//______________________________________________________________________________
	49	AliTagCreator::AliTagCreator() //local mode
	50	{
	51	//==============Default constructor for a AliTagCreator==================
28afeb2e	52	fgridpath = "";
b45e5084	53	fSE = "";
28afeb2e	54	fStorage = 0;
b45e5084	55	}
b45e5084	56
b45e5084	57	//______________________________________________________________________________
	58	AliTagCreator::~AliTagCreator()
	59	{
	60	//================Default destructor for a AliTagCreator=======================
	61	}
	62
28afeb2e	63	//______________________________________________________________________________
	64	void AliTagCreator::SetStorage(Int_t storage)
	65	{
cb1645b7	66	// Sets correctly the storage: 0 for local, 1 for GRID
28afeb2e	67	fStorage = storage;
	68	if(fStorage == 0)
	69	AliInfo(Form("Tags will be stored locally...."));
	70	if(fStorage == 1)
	71	AliInfo(Form("Tags will be stored in the grid...."));
	72	if((fStorage != 0)&&(fStorage != 1))
	73	{
	74	AliInfo(Form("Storage was not properly set!!!"));
	75	abort();
	76	}
	77	}
	78
b45e5084	79
b45e5084	80	//______________________________________________________________________________
cb1645b7	81	Bool_t AliTagCreator::ReadESDCollection(TGridResult *fresult)
b45e5084	82	{
cb1645b7	83	// Reads the entry of the TGridResult and creates the tags
cb1645b7	84	Int_t nEntries = fresult->GetEntries();
b45e5084	85
cb1645b7	86	TString alienUrl;
b45e5084	87	const char *guid;
e16601cf	88	const char *md5;
	89	const char *turl;
	90	Long64_t size = -1;
	91
b45e5084	92	Int_t counter = 0;
cb1645b7	93	for(Int_t i = 0; i < nEntries; i++)
09c52d4a	94	{
cb1645b7	95	alienUrl = fresult->GetKey(i,"turl");
09c52d4a	96	guid = fresult->GetKey(i,"guid");
e16601cf	97	if(fresult->GetKey(i,"size"))
	98	size = atol (fresult->GetKey(i,"size"));
	99	md5 = fresult->GetKey(i,"md5");
	100	turl = fresult->GetKey(i,"turl");
	101	if(md5 && !strlen(guid))
	102	md5 = 0;
	103	if(guid && !strlen(guid))
	104	guid = 0;
	105
cb1645b7	106	TFile *f = TFile::Open(alienUrl,"READ");
e16601cf	107	CreateTag(f,guid,md5,turl,size,counter);
09c52d4a	108	f->Close();
09c52d4a	109	delete f;
b45e5084	110	counter += 1;
09c52d4a	111	}//grid result loop
b45e5084	112
	113	return kTRUE;
	114	}
	115
4fd84456	116	//__________________________________________________________________________
	117	Bool_t AliTagCreator::MergeTags()
	118	{
	119	//Merges the tags and stores the merged tag file
	120	//locally if fStorage=0 or in the grid if fStorage=1
	121	AliInfo(Form("Merging tags....."));
	122	TChain *fgChain = new TChain("T");
	123
	124	if(fStorage == 0) {
	125	const char * tagPattern = "tag";
	126	// Open the working directory
	127	void * dirp = gSystem->OpenDirectory(gSystem->pwd());
	128	const char * name = 0x0;
	129	// Add all files matching pattern to the chain
	130	while((name = gSystem->GetDirEntry(dirp))) {
	131	if (strstr(name,tagPattern))
	132	fgChain->Add(name);
	133	}//directory loop
	134	AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
	135	}//local mode
	136
	137	else if(fStorage == 1) {
	138	TString alienLocation = gGrid->Pwd();
	139	alienLocation += fgridpath.Data();
	140	alienLocation += "/";
	141
	142	TGridResult tagresult = gGrid->Query(alienLocation,"tag.root","","");
	143	Int_t nEntries = tagresult->GetEntries();
	144	for(Int_t i = 0; i < nEntries; i++) {
	145	TString alienUrl = tagresult->GetKey(i,"turl");
	146	fgChain->Add(alienUrl);
	147	}//grid result loop
	148	AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
	149	}//grid mode
	150
	151	AliRunTag *tag = new AliRunTag;
	152	AliEventTag *evTag = new AliEventTag;
	153	fgChain->SetBranchAddress("AliTAG",&tag);
	154
	155	//Defining new tag objects
	156	AliRunTag *newTag = new AliRunTag();
	157	TTree ttag("T","A Tree with event tags");
	158	TBranch * btag = ttag.Branch("AliTAG", &newTag);
	159	btag->SetCompressionLevel(9);
	160	for(Int_t iTagFiles = 0; iTagFiles < fgChain->GetEntries(); iTagFiles++) {
	161	fgChain->GetEntry(iTagFiles);
	162	newTag->SetRunId(tag->GetRunId());
	163	const TClonesArray *tagList = tag->GetEventTags();
	164	for(Int_t j = 0; j < tagList->GetEntries(); j++) {
	165	evTag = (AliEventTag *) tagList->At(j);
	166	newTag->AddEventTag(*evTag);
	167	}
	168	ttag.Fill();
	169	newTag->Clear();
	170	}//tag file loop
	171
	172	TString localFileName = "Run"; localFileName += tag->GetRunId();
	173	localFileName += ".Merged"; localFileName += ".ESD.tag.root";
	174
	175	TString alienFileName = "/alien";
	176	alienFileName += gGrid->Pwd();
	177	alienFileName += fgridpath.Data();
	178	alienFileName += "/";
	179	alienFileName += localFileName;
180	alienFileName += "?se=";
181	alienFileName += fSE.Data();
182
183	TString filename = 0x0;
184
185	if(fStorage == 0) {
186	filename = localFileName.Data();
187	AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
188	}
189	if(fStorage == 1) {
190	filename = alienFileName.Data();
191	AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));
192	}
193
194	TFile* ftag = TFile::Open(filename, "recreate");
195	ftag->cd();
196	ttag.Write();
197	ftag->Close();
198
199	delete tag;
4fd84456	200	delete newTag;
	201
	202	return kTRUE;
	203	}
b45e5084	204
b45e5084	205	//_____________________________________________________________________________
e16601cf	206	void AliTagCreator::CreateTag(TFile* file, const char guid, const char md5, const char *turl, Long64_t size, Int_t Counter)
b45e5084	207	{
2bdb9d38	208	/////////////
	209	//muon code//
	210	////////////
56982dd1	211	Double_t fMUONMASS = 0.105658369;
2bdb9d38	212	//Variables
56982dd1	213	Double_t fX,fY,fZ ;
	214	Double_t fThetaX, fThetaY, fPyz, fChisquare;
	215	Double_t fPxRec, fPyRec, fPzRec, fEnergy;
	216	Int_t fCharge;
	217	TLorentzVector fEPvector;
	218
	219	Float_t fZVertexCut = 10.0;
	220	Float_t fRhoVertexCut = 2.0;
	221
	222	Float_t fLowPtCut = 1.0;
	223	Float_t fHighPtCut = 3.0;
	224	Float_t fVeryHighPtCut = 10.0;
2bdb9d38	225	////////////
	226
	227	Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
	228
cb1645b7	229	// Creates the tags for all the events in a given ESD file
b45e5084	230	Int_t ntrack;
cb1645b7	231	Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
	232	Int_t nPos, nNeg, nNeutr;
	233	Int_t nK0s, nNeutrons, nPi0s, nGamas;
	234	Int_t nCh1GeV, nCh3GeV, nCh10GeV;
	235	Int_t nMu1GeV, nMu3GeV, nMu10GeV;
	236	Int_t nEl1GeV, nEl3GeV, nEl10GeV;
	237	Float_t maxPt = .0, meanPt = .0, totalP = .0;
d1a50cb5	238	Int_t fVertexflag;
d1a50cb5	239	TString fVertexName;
b45e5084	240
b45e5084	241	AliRunTag *tag = new AliRunTag();
b45e5084	242	AliEventTag *evTag = new AliEventTag();
b45e5084	243	TTree ttag("T","A Tree with event tags");
4fd84456	244	TBranch * btag = ttag.Branch("AliTAG", &tag);
b45e5084	245	btag->SetCompressionLevel(9);
	246
	247	AliInfo(Form("Creating the tags......."));
	248
	249	Int_t firstEvent = 0,lastEvent = 0;
	250	TTree t = (TTree) file->Get("esdTree");
	251	TBranch * b = t->GetBranch("ESD");
	252	AliESD *esd = 0;
	253	b->SetAddress(&esd);
	254
	255	tag->SetRunId(esd->GetRunNumber());
	256
cb1645b7	257	Int_t iNumberOfEvents = b->GetEntries();
2bdb9d38	258	for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
	259	ntrack = 0;
	260	nPos = 0;
	261	nNeg = 0;
	262	nNeutr =0;
	263	nK0s = 0;
	264	nNeutrons = 0;
	265	nPi0s = 0;
	266	nGamas = 0;
	267	nProtons = 0;
	268	nKaons = 0;
	269	nPions = 0;
	270	nMuons = 0;
	271	nElectrons = 0;
	272	nCh1GeV = 0;
	273	nCh3GeV = 0;
	274	nCh10GeV = 0;
	275	nMu1GeV = 0;
	276	nMu3GeV = 0;
	277	nMu10GeV = 0;
	278	nEl1GeV = 0;
	279	nEl3GeV = 0;
	280	nEl10GeV = 0;
	281	maxPt = .0;
	282	meanPt = .0;
	283	totalP = .0;
d1a50cb5	284	fVertexflag = 1;
2bdb9d38	285
	286	b->GetEntry(iEventNumber);
	287	const AliESDVertex * vertexIn = esd->GetVertex();
d1a50cb5	288	fVertexName = vertexIn->GetName();
	289	if(fVertexName == "default") fVertexflag = 0;
	290
2bdb9d38	291	for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
	292	AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
	293	UInt_t status = esdTrack->GetStatus();
b45e5084	294
2bdb9d38	295	//select only tracks with ITS refit
	296	if ((status&AliESDtrack::kITSrefit)==0) continue;
	297	//select only tracks with TPC refit
	298	if ((status&AliESDtrack::kTPCrefit)==0) continue;
b45e5084	299
2bdb9d38	300	//select only tracks with the "combined PID"
	301	if ((status&AliESDtrack::kESDpid)==0) continue;
	302	Double_t p[3];
	303	esdTrack->GetPxPyPz(p);
	304	Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
	305	Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
	306	totalP += momentum;
	307	meanPt += fPt;
	308	if(fPt > maxPt) maxPt = fPt;
b45e5084	309
2bdb9d38	310	if(esdTrack->GetSign() > 0) {
2bdb9d38	311	nPos++;
56982dd1	312	if(fPt > fLowPtCut) nCh1GeV++;
	313	if(fPt > fHighPtCut) nCh3GeV++;
	314	if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38	315	}
	316	if(esdTrack->GetSign() < 0) {
	317	nNeg++;
56982dd1	318	if(fPt > fLowPtCut) nCh1GeV++;
	319	if(fPt > fHighPtCut) nCh3GeV++;
	320	if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38	321	}
2bdb9d38	322	if(esdTrack->GetSign() == 0) nNeutr++;
b45e5084	323
2bdb9d38	324	//PID
	325	Double_t prob[5];
	326	esdTrack->GetESDpid(prob);
b45e5084	327
2bdb9d38	328	Double_t rcc = 0.0;
	329	for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
	330	if(rcc == 0.0) continue;
	331	//Bayes' formula
	332	Double_t w[5];
	333	for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
b45e5084	334
2bdb9d38	335	//protons
	336	if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
	337	//kaons
	338	if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
	339	//pions
	340	if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
	341	//electrons
	342	if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
	343	nElectrons++;
56982dd1	344	if(fPt > fLowPtCut) nEl1GeV++;
	345	if(fPt > fHighPtCut) nEl3GeV++;
	346	if(fPt > fVeryHighPtCut) nEl10GeV++;
2bdb9d38	347	}
	348	ntrack++;
	349	}//esd track loop
	350
	351	/////////////
	352	//muon code//
	353	////////////
	354	Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
	355	// loop over all reconstructed tracks (also first track of combination)
	356	for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
	357	AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
	358	if (muonTrack == 0x0) continue;
b45e5084	359
2bdb9d38	360	// Coordinates at vertex
56982dd1	361	fZ = muonTrack->GetZ();
	362	fY = muonTrack->GetBendingCoor();
	363	fX = muonTrack->GetNonBendingCoor();
b45e5084	364
56982dd1	365	fThetaX = muonTrack->GetThetaX();
56982dd1	366	fThetaY = muonTrack->GetThetaY();
b45e5084	367
56982dd1	368	fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
	369	fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
	370	fPxRec = fPzRec * TMath::Tan(fThetaX);
	371	fPyRec = fPzRec * TMath::Tan(fThetaY);
	372	fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
b45e5084	373
2bdb9d38	374	//ChiSquare of the track if needed
56982dd1	375	fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
	376	fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
	377	fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
b45e5084	378
2bdb9d38	379	// total number of muons inside a vertex cut
56982dd1	380	if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fYfY+fXfX)<fRhoVertexCut)) {
2bdb9d38	381	nMuons++;
56982dd1	382	if(fEPvector.Pt() > fLowPtCut) {
2bdb9d38	383	nMu1GeV++;
56982dd1	384	if(fEPvector.Pt() > fHighPtCut) {
2bdb9d38	385	nMu3GeV++;
56982dd1	386	if (fEPvector.Pt() > fVeryHighPtCut) {
2bdb9d38	387	nMu10GeV++;
	388	}
	389	}
	390	}
	391	}
	392	}//muon track loop
	393
	394	// Fill the event tags
	395	if(ntrack != 0)
	396	meanPt = meanPt/ntrack;
	397
	398	evTag->SetEventId(iEventNumber+1);
	399	evTag->SetGUID(guid);
	400	evTag->SetMD5(md5);
	401	evTag->SetTURL(turl);
	402	evTag->SetSize(size);
	403	evTag->SetVertexX(vertexIn->GetXv());
	404	evTag->SetVertexY(vertexIn->GetYv());
	405	evTag->SetVertexZ(vertexIn->GetZv());
d1a50cb5	406	evTag->SetVertexZError(vertexIn->GetZRes());
d1a50cb5	407	evTag->SetVertexFlag(fVertexflag);
2bdb9d38	408
	409	evTag->SetT0VertexZ(esd->GetT0zVertex());
	410
b647652d	411	evTag->SetTrigger(esd->GetTriggerMask());
2bdb9d38	412
32a5cab4	413	evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
32a5cab4	414	evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
2bdb9d38	415	evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
32a5cab4	416	evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
32a5cab4	417	evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
2bdb9d38	418	evTag->SetNumOfParticipants(esd->GetZDCParticipants());
	419
	420
	421	evTag->SetNumOfTracks(esd->GetNumberOfTracks());
	422	evTag->SetNumOfPosTracks(nPos);
	423	evTag->SetNumOfNegTracks(nNeg);
	424	evTag->SetNumOfNeutrTracks(nNeutr);
	425
	426	evTag->SetNumOfV0s(esd->GetNumberOfV0s());
	427	evTag->SetNumOfCascades(esd->GetNumberOfCascades());
	428	evTag->SetNumOfKinks(esd->GetNumberOfKinks());
	429	evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
	430
	431	evTag->SetNumOfProtons(nProtons);
	432	evTag->SetNumOfKaons(nKaons);
	433	evTag->SetNumOfPions(nPions);
	434	evTag->SetNumOfMuons(nMuons);
	435	evTag->SetNumOfElectrons(nElectrons);
	436	evTag->SetNumOfPhotons(nGamas);
	437	evTag->SetNumOfPi0s(nPi0s);
	438	evTag->SetNumOfNeutrons(nNeutrons);
	439	evTag->SetNumOfKaon0s(nK0s);
	440
	441	evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
	442	evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
	443	evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
	444	evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
	445	evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
	446	evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
	447	evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
	448	evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
	449	evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
	450
85c60a8e	451	evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
85c60a8e	452	evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
2bdb9d38	453
	454	evTag->SetTotalMomentum(totalP);
	455	evTag->SetMeanPt(meanPt);
	456	evTag->SetMaxPt(maxPt);
	457
	458	tag->AddEventTag(*evTag);
	459	}//event loop
cb1645b7	460	lastEvent = iNumberOfEvents;
b45e5084	461
b45e5084	462	t->Delete("");
2bdb9d38	463
b45e5084	464	ttag.Fill();
	465	tag->Clear();
	466
cb1645b7	467	TString localFileName = "Run"; localFileName += tag->GetRunId();
	468	localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
	469	localFileName += ".ESD.tag.root";
b45e5084	470
cb1645b7	471	TString alienLocation = "/alien";
	472	alienLocation += gGrid->Pwd();
	473	alienLocation += fgridpath.Data();
	474	alienLocation += "/";
	475	alienLocation += localFileName;
	476	alienLocation += "?se=";
	477	alienLocation += fSE.Data();
28afeb2e	478
cb1645b7	479	TString fileName;
28afeb2e	480
	481	if(fStorage == 0)
	482	{
cb1645b7	483	fileName = localFileName.Data();
cb1645b7	484	AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
28afeb2e	485	}
	486	if(fStorage == 1)
	487	{
cb1645b7	488	fileName = alienLocation.Data();
cb1645b7	489	AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
28afeb2e	490	}
b45e5084	491
cb1645b7	492	TFile* ftag = TFile::Open(fileName, "recreate");
b45e5084	493	ftag->cd();
	494	ttag.Write();
	495	ftag->Close();
	496
	497	delete ftag;
	498	delete esd;
	499
	500	delete tag;
b45e5084	501	}
b45e5084	502