[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliJetEmbeddingFromAODTask.cxx

// $Id: AliJetEmbeddingFromAODTask.cxx $
//
// Jet embedding from AOD task.
//
// Author: S.Aiola, C.Loizides

#include "AliJetEmbeddingFromAODTask.h"

#include <TFile.h>
#include <TTree.h>
#include <TClonesArray.h>
#include <TObjArray.h>
#include <TObjString.h>
#include <TGrid.h>
#include <TH2C.h>
#include <TList.h>
#include <TStreamerInfo.h>
#include <TRandom.h>

#include "AliVEvent.h"
#include "AliAODTrack.h"
#include "AliESDtrack.h"
#include "AliPicoTrack.h"
#include "AliVTrack.h"
#include "AliVCluster.h"
#include "AliVCaloCells.h"
#include "AliAODMCParticle.h"
#include "AliCentrality.h"
#include "AliVHeader.h"
#include "AliVVertex.h"
#include "AliAODHeader.h"
#include "AliLog.h"
#include "AliInputEventHandler.h"

ClassImp(AliJetEmbeddingFromAODTask)

//________________________________________________________________________
AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask() : 
  AliJetModelBaseTask("AliJetEmbeddingFromAODTask"),
  fFileList(0),
  fRandomAccess(kFALSE),
  fAODTreeName(),
  fAODHeaderName(),
  fAODVertexName(),
  fAODTrackName(),
  fAODClusName(),
  fAODCellsName(),
  fAODMCParticlesName(),
  fMinCentrality(0),
  fMaxCentrality(10),
  fTriggerMask(AliVEvent::kAny),
  fZVertexCut(10),
  fIncludeNoITS(kTRUE),
  fTotalFiles(2050),
  fAttempts(5),
  fEsdTreeMode(kFALSE),
  fCurrentFileID(0),
  fCurrentAODFileID(0),
  fCurrentAODFile(0),
  fPicoTrackVersion(0),
  fCurrentAODTree(0),
  fAODHeader(0),
  fAODVertex(0),
  fAODTracks(0),
  fAODClusters(0),
  fAODCaloCells(0),
  fAODMCParticles(0),
  fCurrentAODEntry(-1),
  fHistFileMatching(0),
  fHistAODFileError(0),
  fHistNotEmbedded(0)
{
  // Default constructor.
  SetSuffix("AODEmbedding");
  SetMarkMC(0);
  fAODfilterBits[0] = -1;
  fAODfilterBits[1] = -1;
}

//________________________________________________________________________
AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask(const char *name, Bool_t drawqa) : 
  AliJetModelBaseTask(name, drawqa),
  fFileList(0),
  fRandomAccess(kFALSE),
  fAODTreeName("aodTree"),
  fAODHeaderName("header"),
  fAODVertexName("vertices"),
  fAODTrackName("tracks"),
  fAODClusName(""),
  fAODCellsName("emcalCells"),
  fAODMCParticlesName(AliAODMCParticle::StdBranchName()),
  fMinCentrality(0),
  fMaxCentrality(10),
  fTriggerMask(AliVEvent::kAny),
  fZVertexCut(10),
  fIncludeNoITS(kTRUE),
  fTotalFiles(2050),
  fAttempts(5),
  fEsdTreeMode(kFALSE),
  fCurrentFileID(0),
  fCurrentAODFileID(0),
  fCurrentAODFile(0),
  fPicoTrackVersion(0),
  fCurrentAODTree(0),
  fAODHeader(0),
  fAODVertex(0),
  fAODTracks(0),
  fAODClusters(0),
  fAODCaloCells(0),  
  fAODMCParticles(0),
  fCurrentAODEntry(-1),
  fHistFileMatching(0),
  fHistAODFileError(0),
  fHistNotEmbedded(0)
{
  // Standard constructor.
  SetSuffix("AODEmbedding");
  SetMarkMC(0);
  fAODfilterBits[0] = -1;
  fAODfilterBits[1] = -1;
}

//________________________________________________________________________
AliJetEmbeddingFromAODTask::~AliJetEmbeddingFromAODTask()
{
  // Destructor

  if (fCurrentAODFile) {
    fCurrentAODFile->Close();
    delete fCurrentAODFile;
  }
}

//________________________________________________________________________
void AliJetEmbeddingFromAODTask::UserCreateOutputObjects()
{
  if (!fQAhistos)
    return;

  AliJetModelBaseTask::UserCreateOutputObjects();
  
  fHistFileMatching = new TH2C("fHistFileMatching", "fHistFileMatching", fTotalFiles, -0.5, fTotalFiles - 0.5, fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
  fHistFileMatching->GetXaxis()->SetTitle("File no. (PYTHIA)");
  fHistFileMatching->GetYaxis()->SetTitle("File no. (Embedded AOD)");
  fHistFileMatching->GetZaxis()->SetTitle("counts");
  fOutput->Add(fHistFileMatching);

  fHistAODFileError = new TH1C("fHistAODFileError", "fHistAODFileError", fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
  fHistAODFileError->GetXaxis()->SetTitle("File no. (Embedded AOD)");
  fHistAODFileError->GetYaxis()->SetTitle("counts");
  fOutput->Add(fHistAODFileError);

  fHistNotEmbedded = new TH1C("fHistNotEmbedded", "fHistNotEmbedded", fTotalFiles, -0.5, fTotalFiles - 0.5);
  fHistNotEmbedded->GetXaxis()->SetTitle("File no. (PYTHIA)");
  fHistNotEmbedded->GetYaxis()->SetTitle("counts");
  fOutput->Add(fHistNotEmbedded);

  PostData(1, fOutput);
}

//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::UserNotify()
{
  TString path(fInputHandler->GetTree()->GetCurrentFile()->GetPath());
  if (path.EndsWith("/"))
    path.Remove(path.Length()-1);
  path.Remove(path.Last('/'));
  path.Remove(0, path.Last('/')+1);
  if (!path.IsDec()) {
    AliWarning(Form("Could not extract file number from path %s", path.Data()));
    return kTRUE;
  }
  fCurrentFileID = path.Atoi();
  if (!fRandomAccess) {
    fCurrentAODFileID = fFileList->GetEntriesFast() * fCurrentFileID / fTotalFiles-1;
    AliInfo(Form("Start embedding from file ID %d", fCurrentAODFileID));
  }
  return kTRUE;
}

//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::ExecOnce() 
{
  if (fAODTreeName.Contains("aod"))
    fEsdTreeMode = kFALSE;
  else
    fEsdTreeMode = kTRUE;

  return AliJetModelBaseTask::ExecOnce();
}

//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::OpenNextFile() 
{
  if (fCurrentAODFile) {
    fCurrentAODFile->Close();
    delete fCurrentAODFile;
    fCurrentAODFile = 0;
  }

  Int_t i = 0;
  TString fileName;

  do {
    if (i>0) {
      AliDebug(3,Form("Failed to open file %s...", fileName.Data()));
      if (fHistAODFileError)
	fHistAODFileError->Fill(fCurrentAODFileID);
    }

    fileName = GetNextFileName();

    if (fileName.IsNull())
      break;
    
    if (fileName.BeginsWith("alien://") && !gGrid) {
      AliInfo("Trying to connect to AliEn ...");
      TGrid::Connect("alien://");
    }

    AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
    fCurrentAODFile = TFile::Open(fileName);

    i++;

  } while ((!fCurrentAODFile || fCurrentAODFile->IsZombie()) && i < fAttempts);

  if (!fCurrentAODFile || fCurrentAODFile->IsZombie())
    return kFALSE;

  const TList *clist = fCurrentAODFile->GetStreamerInfoCache();
  if(clist) {
    TStreamerInfo *cinfo = static_cast<TStreamerInfo*>(clist->FindObject("AliPicoTrack"));
    if(cinfo) 
      fPicoTrackVersion = cinfo->GetClassVersion();
    else
      fPicoTrackVersion = 0;
  }

  fCurrentAODTree = static_cast<TTree*>(fCurrentAODFile->Get(fAODTreeName));
  if (!fCurrentAODTree)
    return kFALSE;

  if (!fAODHeaderName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODHeaderName, &fAODHeader);
  
  if (!fAODVertexName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODVertexName, &fAODVertex);
      
  if (!fAODTrackName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODTrackName, &fAODTracks);
  
  if (!fAODClusName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODClusName, &fAODClusters);
  
  if (!fAODCellsName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODCellsName, &fAODCaloCells);
  
  if (!fAODMCParticlesName.IsNull()) 
    fCurrentAODTree->SetBranchAddress(fAODMCParticlesName, &fAODMCParticles);
  
  if (fRandomAccess)
    fCurrentAODEntry = TMath::Nint(gRandom->Rndm()*fCurrentAODTree->GetEntries());
  else
    fCurrentAODEntry = 0;

  AliDebug(2,Form("Will start embedding from entry %d", fCurrentAODEntry));
  
  if (fHistFileMatching)
    fHistFileMatching->Fill(fCurrentFileID, fCurrentAODFileID-1);
  
  return kTRUE;
}

//________________________________________________________________________
TString AliJetEmbeddingFromAODTask::GetNextFileName()
{
    if (fRandomAccess) 
      fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*fFileList->GetEntriesFast());
    else
      fCurrentAODFileID++;

    if (fCurrentAODFileID >= fFileList->GetEntriesFast())
      return "";
    
    TObjString *objFileName = static_cast<TObjString*>(fFileList->At(fCurrentAODFileID));
    return objFileName->GetString();
}

//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::GetNextEntry() 
{
  Int_t attempts = 0;

  do {
    if (!fCurrentAODFile || !fCurrentAODTree || fCurrentAODEntry >= fCurrentAODTree->GetEntries()) {
      if (!OpenNextFile())
	return kFALSE;
    }
    
    fCurrentAODTree->GetEntry(fCurrentAODEntry);
    fCurrentAODEntry++;
    attempts++;

    if (attempts == 1000) 
      AliWarning("After 1000 attempts no event has been accepted by the event selection (trigger, centrality...)!");

  } while (!IsAODEventSelected() && fCurrentAODTree);

  return (fCurrentAODTree!=0);
}

//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::IsAODEventSelected()
{
  if (!fEsdTreeMode && fAODHeader) {
    AliAODHeader *aodHeader = static_cast<AliAODHeader*>(fAODHeader);

    if (fTriggerMask != AliVEvent::kAny) {
      UInt_t offlineTrigger = aodHeader->GetOfflineTrigger();
      
      if ((offlineTrigger & fTriggerMask) == 0) {
	AliDebug(2,Form("Event rejected due to physics selection. Event trigger mask: %d, trigger mask selection: %d.", offlineTrigger, fTriggerMask));
	return kFALSE;
      }
    }
    
    if (fMinCentrality >= 0) {
      AliCentrality *cent = aodHeader->GetCentralityP();
      Float_t centVal = cent->GetCentralityPercentile("V0M");
      if (centVal < fMinCentrality || centVal >= fMaxCentrality) {
	AliDebug(2,Form("Event rejected due to centrality selection. Event centrality: %f, centrality range selection: %f to %f", centVal, fMinCentrality, fMaxCentrality));
	return kFALSE;
      }
    }
  }

  if (fAODVertex) {
    Double_t vert[3]={0};
    ((AliVVertex*)fAODVertex->At(0))->GetXYZ(vert);
    if (TMath::Abs(vert[2]) > fZVertexCut) {
      AliDebug(2,Form("Event rejected due to Z vertex selection. Event Z vertex: %f, Z vertex cut: %f", vert[2], fZVertexCut));
      return kFALSE;
    }
  }

  return kTRUE;
}

//________________________________________________________________________
void AliJetEmbeddingFromAODTask::Run() 
{
  if (!GetNextEntry()) {
    if (fHistNotEmbedded)
      fHistNotEmbedded->Fill(fCurrentFileID);
    AliError("Unable to get the AOD event to embed. Nothing will be embedded.");
    return;
  }

  if (fOutMCParticles) {

    if (fCopyArray && fMCParticles)
      CopyMCParticles();

    if (fAODMCParticles) {
      for (Int_t i = 0; i < fAODMCParticles->GetEntriesFast(); i++) {
	AliAODMCParticle *part = static_cast<AliAODMCParticle*>(fAODMCParticles->At(i));
	if (!part) {
	  AliError(Form("Could not find MC particle %d in branch %s of tree %s!", i, fAODMCParticlesName.Data(), fAODTreeName.Data()));
	  continue;
	}
	
	if (!part->IsPhysicalPrimary()) 
	  continue;

	AliDebug(3, Form("Embedding MC particle with pT = %f, eta = %f, phi = %f", part->Pt(), part->Eta(), part->Phi()));
	AddMCParticle(part, i);
      }
    }
  }

  if (fOutTracks) {

    if (fCopyArray && fTracks)
      CopyTracks();

    if (fAODTracks) {
      AliDebug(2, Form("%d tracks will be processed for embedding.", fAODTracks->GetEntriesFast()));
      for (Int_t i = 0; i < fAODTracks->GetEntriesFast(); i++) {
	AliVTrack *track = static_cast<AliVTrack*>(fAODTracks->At(i));
	if (!track) {
	  AliError(Form("Could not find track %d in branch %s of tree %s!", i, fAODTrackName.Data(), fAODTreeName.Data()));
	  continue;
	}
	
	Int_t type = 0;
	Bool_t isEmc = kFALSE;
	if (!fEsdTreeMode) {
	  AliAODTrack *aodtrack = static_cast<AliAODTrack*>(track);
	  if (aodtrack->TestFilterBit(fAODfilterBits[0])) {
	    type = 0;
	  }
	  else if (aodtrack->TestFilterBit(fAODfilterBits[1])) {
	    if ((aodtrack->GetStatus()&AliESDtrack::kITSrefit)==0) {
	      if (fIncludeNoITS)
		type = 2;
	      else
		continue;
	    }
	    else {
	      type = 1;
	    }
	  }
	  else { /*not a good track*/
	    continue;
	  }

	  if (TMath::Abs(aodtrack->GetTrackEtaOnEMCal()) < 0.75 && 
	      aodtrack->GetTrackPhiOnEMCal() > 70 * TMath::DegToRad() && 
	      aodtrack->GetTrackPhiOnEMCal() < 190 * TMath::DegToRad())
	    isEmc = kTRUE;
	}
	else if (fPicoTrackVersion > 0) { /*not AOD mode, let's see if it is PicoTrack*/
	  AliPicoTrack *ptrack = dynamic_cast<AliPicoTrack*>(track);
	  if (ptrack) {
	    if (fPicoTrackVersion >= 3)
	      type = ptrack->GetTrackType();
	    else
	      type = ptrack->GetLabel();
	    isEmc = ptrack->IsEMCAL();
	  }
	}
	
	Int_t label = 0;

	if (fOutMCParticles) {
	  label = track->GetLabel();

	  if (label == 0)
	    label = 99999;
	}

	AliDebug(3, Form("Embedding track with pT = %f, eta = %f, phi = %f, label = %d", track->Pt(), track->Eta(), track->Phi(), label));
	AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc, label);
      }
    }
  }

  if (fOutClusters) {

    if (fCopyArray && fClusters)
      CopyClusters();

    if (fAODClusters) {
      for (Int_t i = 0; i < fAODClusters->GetEntriesFast(); i++) {
	AliVCluster *clus = static_cast<AliVCluster*>(fAODClusters->At(i));
	if (!clus) {
	  AliError(Form("Could not find cluster %d in branch %s of tree %s!", i, fAODClusName.Data(), fAODTreeName.Data()));
	  continue;
	}
	TLorentzVector vect;
	Double_t vert[3] = {0,0,0};
	clus->GetMomentum(vect,vert);
	AddCluster(clus->E(), vect.Eta(), vect.Phi(), clus->GetLabel());
      }
    }
  }

  if (fOutCaloCells) {

    Int_t totalCells = 0;

    if (fCaloCells)
      totalCells += fCaloCells->GetNumberOfCells();
    if (fAODCaloCells)
      totalCells += fAODCaloCells->GetNumberOfCells();

    SetNumberOfOutCells(totalCells);

    if (fCaloCells)
      CopyCells();

    if (fAODCaloCells) {
      for (Short_t i = 0; i < fAODCaloCells->GetNumberOfCells(); i++) {
	Int_t mclabel = 0;
	Double_t efrac = 0.;
	Double_t time = -1;
	Short_t cellNum = -1;
	Double_t amp = -1;
	
	fAODCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
	AliDebug(3,Form("Adding cell with amplitude %f, absolute ID %d, time %f", amp, cellNum, time));
	AddCell(amp, cellNum, time, mclabel);
      }
    }

    AliDebug(2,Form("Added cells = %d, total cells = %d", fAddedCells, totalCells));
  }
}
Commit	Line	Data
b16bb001	1	// $Id: AliJetEmbeddingFromAODTask.cxx $
	2	//
	3	// Jet embedding from AOD task.
	4	//
	5	// Author: S.Aiola, C.Loizides
	6
	7	#include "AliJetEmbeddingFromAODTask.h"
	8
b16bb001	9	#include <TFile.h>
	10	#include <TTree.h>
	11	#include <TClonesArray.h>
	12	#include <TObjArray.h>
	13	#include <TObjString.h>
	14	#include <TGrid.h>
	15	#include <TH2C.h>
cd6431de	16	#include <TList.h>
cd6431de	17	#include <TStreamerInfo.h>
6a20534a	18	#include <TRandom.h>
b16bb001	19
	20	#include "AliVEvent.h"
	21	#include "AliAODTrack.h"
	22	#include "AliESDtrack.h"
	23	#include "AliPicoTrack.h"
	24	#include "AliVTrack.h"
	25	#include "AliVCluster.h"
	26	#include "AliVCaloCells.h"
787a3c4f	27	#include "AliAODMCParticle.h"
b16bb001	28	#include "AliCentrality.h"
	29	#include "AliVHeader.h"
	30	#include "AliVVertex.h"
	31	#include "AliAODHeader.h"
	32	#include "AliLog.h"
	33	#include "AliInputEventHandler.h"
	34
	35	ClassImp(AliJetEmbeddingFromAODTask)
	36
	37	//________________________________________________________________________
	38	AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask() :
	39	AliJetModelBaseTask("AliJetEmbeddingFromAODTask"),
	40	fFileList(0),
6a20534a	41	fRandomAccess(kFALSE),
b16bb001	42	fAODTreeName(),
	43	fAODHeaderName(),
	44	fAODVertexName(),
	45	fAODTrackName(),
	46	fAODClusName(),
	47	fAODCellsName(),
787a3c4f	48	fAODMCParticlesName(),
b16bb001	49	fMinCentrality(0),
	50	fMaxCentrality(10),
	51	fTriggerMask(AliVEvent::kAny),
	52	fZVertexCut(10),
	53	fIncludeNoITS(kTRUE),
6a20534a	54	fTotalFiles(2050),
6a20534a	55	fAttempts(5),
b16bb001	56	fEsdTreeMode(kFALSE),
b16bb001	57	fCurrentFileID(0),
6a20534a	58	fCurrentAODFileID(0),
b16bb001	59	fCurrentAODFile(0),
cd6431de	60	fPicoTrackVersion(0),
6a20534a	61	fCurrentAODTree(0),
b16bb001	62	fAODHeader(0),
	63	fAODVertex(0),
	64	fAODTracks(0),
	65	fAODClusters(0),
	66	fAODCaloCells(0),
787a3c4f	67	fAODMCParticles(0),
6a20534a	68	fCurrentAODEntry(-1),
	69	fHistFileMatching(0),
	70	fHistAODFileError(0),
	71	fHistNotEmbedded(0)
b16bb001	72	{
	73	// Default constructor.
	74	SetSuffix("AODEmbedding");
7f76e479	75	SetMarkMC(0);
b16bb001	76	fAODfilterBits[0] = -1;
	77	fAODfilterBits[1] = -1;
	78	}
	79
	80	//________________________________________________________________________
	81	AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask(const char *name, Bool_t drawqa) :
	82	AliJetModelBaseTask(name, drawqa),
	83	fFileList(0),
6a20534a	84	fRandomAccess(kFALSE),
b16bb001	85	fAODTreeName("aodTree"),
	86	fAODHeaderName("header"),
	87	fAODVertexName("vertices"),
	88	fAODTrackName("tracks"),
787a3c4f	89	fAODClusName(""),
b16bb001	90	fAODCellsName("emcalCells"),
787a3c4f	91	fAODMCParticlesName(AliAODMCParticle::StdBranchName()),
b16bb001	92	fMinCentrality(0),
	93	fMaxCentrality(10),
	94	fTriggerMask(AliVEvent::kAny),
	95	fZVertexCut(10),
	96	fIncludeNoITS(kTRUE),
6a20534a	97	fTotalFiles(2050),
6a20534a	98	fAttempts(5),
b16bb001	99	fEsdTreeMode(kFALSE),
b16bb001	100	fCurrentFileID(0),
6a20534a	101	fCurrentAODFileID(0),
b16bb001	102	fCurrentAODFile(0),
cd6431de	103	fPicoTrackVersion(0),
6a20534a	104	fCurrentAODTree(0),
b16bb001	105	fAODHeader(0),
	106	fAODVertex(0),
	107	fAODTracks(0),
	108	fAODClusters(0),
787a3c4f	109	fAODCaloCells(0),
787a3c4f	110	fAODMCParticles(0),
6a20534a	111	fCurrentAODEntry(-1),
	112	fHistFileMatching(0),
	113	fHistAODFileError(0),
	114	fHistNotEmbedded(0)
b16bb001	115	{
	116	// Standard constructor.
	117	SetSuffix("AODEmbedding");
7f76e479	118	SetMarkMC(0);
b16bb001	119	fAODfilterBits[0] = -1;
	120	fAODfilterBits[1] = -1;
	121	}
	122
	123	//________________________________________________________________________
	124	AliJetEmbeddingFromAODTask::~AliJetEmbeddingFromAODTask()
	125	{
	126	// Destructor
	127
	128	if (fCurrentAODFile) {
	129	fCurrentAODFile->Close();
	130	delete fCurrentAODFile;
	131	}
	132	}
	133
	134	//________________________________________________________________________
	135	void AliJetEmbeddingFromAODTask::UserCreateOutputObjects()
	136	{
	137	if (!fQAhistos)
	138	return;
	139
	140	AliJetModelBaseTask::UserCreateOutputObjects();
	141
6a20534a	142	fHistFileMatching = new TH2C("fHistFileMatching", "fHistFileMatching", fTotalFiles, -0.5, fTotalFiles - 0.5, fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
	143	fHistFileMatching->GetXaxis()->SetTitle("File no. (PYTHIA)");
	144	fHistFileMatching->GetYaxis()->SetTitle("File no. (Embedded AOD)");
	145	fHistFileMatching->GetZaxis()->SetTitle("counts");
	146	fOutput->Add(fHistFileMatching);
	147
	148	fHistAODFileError = new TH1C("fHistAODFileError", "fHistAODFileError", fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
	149	fHistAODFileError->GetXaxis()->SetTitle("File no. (Embedded AOD)");
	150	fHistAODFileError->GetYaxis()->SetTitle("counts");
	151	fOutput->Add(fHistAODFileError);
	152
	153	fHistNotEmbedded = new TH1C("fHistNotEmbedded", "fHistNotEmbedded", fTotalFiles, -0.5, fTotalFiles - 0.5);
	154	fHistNotEmbedded->GetXaxis()->SetTitle("File no. (PYTHIA)");
	155	fHistNotEmbedded->GetYaxis()->SetTitle("counts");
	156	fOutput->Add(fHistNotEmbedded);
b16bb001	157
	158	PostData(1, fOutput);
	159	}
	160
	161	//________________________________________________________________________
	162	Bool_t AliJetEmbeddingFromAODTask::UserNotify()
	163	{
	164	TString path(fInputHandler->GetTree()->GetCurrentFile()->GetPath());
	165	if (path.EndsWith("/"))
	166	path.Remove(path.Length()-1);
	167	path.Remove(path.Last('/'));
	168	path.Remove(0, path.Last('/')+1);
	169	if (!path.IsDec()) {
	170	AliWarning(Form("Could not extract file number from path %s", path.Data()));
	171	return kTRUE;
	172	}
	173	fCurrentFileID = path.Atoi();
6a20534a	174	if (!fRandomAccess) {
	175	fCurrentAODFileID = fFileList->GetEntriesFast() * fCurrentFileID / fTotalFiles-1;
	176	AliInfo(Form("Start embedding from file ID %d", fCurrentAODFileID));
	177	}
b16bb001	178	return kTRUE;
	179	}
	180
	181	//________________________________________________________________________
	182	Bool_t AliJetEmbeddingFromAODTask::ExecOnce()
	183	{
	184	if (fAODTreeName.Contains("aod"))
	185	fEsdTreeMode = kFALSE;
	186	else
	187	fEsdTreeMode = kTRUE;
	188
	189	return AliJetModelBaseTask::ExecOnce();
	190	}
	191
	192	//________________________________________________________________________
	193	Bool_t AliJetEmbeddingFromAODTask::OpenNextFile()
	194	{
	195	if (fCurrentAODFile) {
	196	fCurrentAODFile->Close();
	197	delete fCurrentAODFile;
	198	fCurrentAODFile = 0;
	199	}
	200
6a20534a	201	Int_t i = 0;
6a20534a	202	TString fileName;
b16bb001	203
6a20534a	204	do {
6a20534a	205	if (i>0) {
7030f36f	206	AliDebug(3,Form("Failed to open file %s...", fileName.Data()));
7030f36f	207	if (fHistAODFileError)
6a20534a	208	fHistAODFileError->Fill(fCurrentAODFileID);
	209	}
	210
	211	fileName = GetNextFileName();
	212
	213	if (fileName.IsNull())
	214	break;
	215
	216	if (fileName.BeginsWith("alien://") && !gGrid) {
b16bb001	217	AliInfo("Trying to connect to AliEn ...");
b16bb001	218	TGrid::Connect("alien://");
6a20534a	219	}
6a20534a	220
7030f36f	221	AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
6a20534a	222	fCurrentAODFile = TFile::Open(fileName);
b16bb001	223
7030f36f	224	i++;
7030f36f	225
6a20534a	226	} while ((!fCurrentAODFile \|\| fCurrentAODFile->IsZombie()) && i < fAttempts);
b16bb001	227
6a20534a	228	if (!fCurrentAODFile \|\| fCurrentAODFile->IsZombie())
b16bb001	229	return kFALSE;
b16bb001	230
cd6431de	231	const TList *clist = fCurrentAODFile->GetStreamerInfoCache();
	232	if(clist) {
	233	TStreamerInfo cinfo = static_cast<TStreamerInfo>(clist->FindObject("AliPicoTrack"));
	234	if(cinfo)
	235	fPicoTrackVersion = cinfo->GetClassVersion();
	236	else
	237	fPicoTrackVersion = 0;
	238	}
	239
6a20534a	240	fCurrentAODTree = static_cast<TTree*>(fCurrentAODFile->Get(fAODTreeName));
	241	if (!fCurrentAODTree)
	242	return kFALSE;
	243
	244	if (!fAODHeaderName.IsNull())
	245	fCurrentAODTree->SetBranchAddress(fAODHeaderName, &fAODHeader);
	246
	247	if (!fAODVertexName.IsNull())
	248	fCurrentAODTree->SetBranchAddress(fAODVertexName, &fAODVertex);
	249
	250	if (!fAODTrackName.IsNull())
	251	fCurrentAODTree->SetBranchAddress(fAODTrackName, &fAODTracks);
	252
	253	if (!fAODClusName.IsNull())
	254	fCurrentAODTree->SetBranchAddress(fAODClusName, &fAODClusters);
	255
	256	if (!fAODCellsName.IsNull())
	257	fCurrentAODTree->SetBranchAddress(fAODCellsName, &fAODCaloCells);
	258
	259	if (!fAODMCParticlesName.IsNull())
	260	fCurrentAODTree->SetBranchAddress(fAODMCParticlesName, &fAODMCParticles);
	261
	262	if (fRandomAccess)
	263	fCurrentAODEntry = TMath::Nint(gRandom->Rndm()*fCurrentAODTree->GetEntries());
	264	else
	265	fCurrentAODEntry = 0;
7030f36f	266
7030f36f	267	AliDebug(2,Form("Will start embedding from entry %d", fCurrentAODEntry));
6a20534a	268
7030f36f	269	if (fHistFileMatching)
6a20534a	270	fHistFileMatching->Fill(fCurrentFileID, fCurrentAODFileID-1);
b16bb001	271
b16bb001	272	return kTRUE;
	273	}
	274
	275	//________________________________________________________________________
6a20534a	276	TString AliJetEmbeddingFromAODTask::GetNextFileName()
b16bb001	277	{
6a20534a	278	if (fRandomAccess)
	279	fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*fFileList->GetEntriesFast());
	280	else
	281	fCurrentAODFileID++;
b16bb001	282
6a20534a	283	if (fCurrentAODFileID >= fFileList->GetEntriesFast())
	284	return "";
	285
	286	TObjString objFileName = static_cast<TObjString>(fFileList->At(fCurrentAODFileID));
	287	return objFileName->GetString();
	288	}
	289
	290	//________________________________________________________________________
	291	Bool_t AliJetEmbeddingFromAODTask::GetNextEntry()
	292	{
b16bb001	293	Int_t attempts = 0;
	294
	295	do {
6a20534a	296	if (!fCurrentAODFile \|\| !fCurrentAODTree \|\| fCurrentAODEntry >= fCurrentAODTree->GetEntries()) {
b16bb001	297	if (!OpenNextFile())
b16bb001	298	return kFALSE;
b16bb001	299	}
b16bb001	300
6a20534a	301	fCurrentAODTree->GetEntry(fCurrentAODEntry);
6a20534a	302	fCurrentAODEntry++;
b16bb001	303	attempts++;
	304
	305	if (attempts == 1000)
	306	AliWarning("After 1000 attempts no event has been accepted by the event selection (trigger, centrality...)!");
	307
6a20534a	308	} while (!IsAODEventSelected() && fCurrentAODTree);
b16bb001	309
6a20534a	310	return (fCurrentAODTree!=0);
b16bb001	311	}
	312
	313	//________________________________________________________________________
	314	Bool_t AliJetEmbeddingFromAODTask::IsAODEventSelected()
	315	{
	316	if (!fEsdTreeMode && fAODHeader) {
	317	AliAODHeader aodHeader = static_cast<AliAODHeader>(fAODHeader);
	318
6a20534a	319	if (fTriggerMask != AliVEvent::kAny) {
	320	UInt_t offlineTrigger = aodHeader->GetOfflineTrigger();
	321
	322	if ((offlineTrigger & fTriggerMask) == 0) {
	323	AliDebug(2,Form("Event rejected due to physics selection. Event trigger mask: %d, trigger mask selection: %d.", offlineTrigger, fTriggerMask));
	324	return kFALSE;
	325	}
b16bb001	326	}
b16bb001	327
6a20534a	328	if (fMinCentrality >= 0) {
	329	AliCentrality *cent = aodHeader->GetCentralityP();
	330	Float_t centVal = cent->GetCentralityPercentile("V0M");
	331	if (centVal < fMinCentrality \|\| centVal >= fMaxCentrality) {
	332	AliDebug(2,Form("Event rejected due to centrality selection. Event centrality: %f, centrality range selection: %f to %f", centVal, fMinCentrality, fMaxCentrality));
	333	return kFALSE;
	334	}
b16bb001	335	}
	336	}
	337
	338	if (fAODVertex) {
	339	Double_t vert[3]={0};
	340	((AliVVertex*)fAODVertex->At(0))->GetXYZ(vert);
	341	if (TMath::Abs(vert[2]) > fZVertexCut) {
	342	AliDebug(2,Form("Event rejected due to Z vertex selection. Event Z vertex: %f, Z vertex cut: %f", vert[2], fZVertexCut));
	343	return kFALSE;
	344	}
	345	}
	346
	347	return kTRUE;
	348	}
	349
	350	//________________________________________________________________________
	351	void AliJetEmbeddingFromAODTask::Run()
	352	{
	353	if (!GetNextEntry()) {
7030f36f	354	if (fHistNotEmbedded)
7030f36f	355	fHistNotEmbedded->Fill(fCurrentFileID);
6a20534a	356	AliError("Unable to get the AOD event to embed. Nothing will be embedded.");
b16bb001	357	return;
	358	}
	359
787a3c4f	360	if (fOutMCParticles) {
	361
	362	if (fCopyArray && fMCParticles)
	363	CopyMCParticles();
	364
	365	if (fAODMCParticles) {
	366	for (Int_t i = 0; i < fAODMCParticles->GetEntriesFast(); i++) {
	367	AliAODMCParticle part = static_cast<AliAODMCParticle>(fAODMCParticles->At(i));
	368	if (!part) {
	369	AliError(Form("Could not find MC particle %d in branch %s of tree %s!", i, fAODMCParticlesName.Data(), fAODTreeName.Data()));
	370	continue;
	371	}
	372
	373	if (!part->IsPhysicalPrimary())
	374	continue;
	375
ca5c29fa	376	AliDebug(3, Form("Embedding MC particle with pT = %f, eta = %f, phi = %f", part->Pt(), part->Eta(), part->Phi()));
787a3c4f	377	AddMCParticle(part, i);
	378	}
	379	}
	380	}
	381
	382	if (fOutTracks) {
b16bb001	383
787a3c4f	384	if (fCopyArray && fTracks)
b16bb001	385	CopyTracks();
	386
	387	if (fAODTracks) {
	388	AliDebug(2, Form("%d tracks will be processed for embedding.", fAODTracks->GetEntriesFast()));
	389	for (Int_t i = 0; i < fAODTracks->GetEntriesFast(); i++) {
	390	AliVTrack track = static_cast<AliVTrack>(fAODTracks->At(i));
	391	if (!track) {
	392	AliError(Form("Could not find track %d in branch %s of tree %s!", i, fAODTrackName.Data(), fAODTreeName.Data()));
	393	continue;
	394	}
	395
	396	Int_t type = 0;
	397	Bool_t isEmc = kFALSE;
	398	if (!fEsdTreeMode) {
	399	AliAODTrack aodtrack = static_cast<AliAODTrack>(track);
	400	if (aodtrack->TestFilterBit(fAODfilterBits[0])) {
	401	type = 0;
	402	}
	403	else if (aodtrack->TestFilterBit(fAODfilterBits[1])) {
	404	if ((aodtrack->GetStatus()&AliESDtrack::kITSrefit)==0) {
	405	if (fIncludeNoITS)
	406	type = 2;
	407	else
	408	continue;
	409	}
	410	else {
	411	type = 1;
	412	}
	413	}
	414	else { /not a good track/
	415	continue;
	416	}
	417
	418	if (TMath::Abs(aodtrack->GetTrackEtaOnEMCal()) < 0.75 &&
	419	aodtrack->GetTrackPhiOnEMCal() > 70 * TMath::DegToRad() &&
	420	aodtrack->GetTrackPhiOnEMCal() < 190 * TMath::DegToRad())
	421	isEmc = kTRUE;
	422	}
cd6431de	423	else if (fPicoTrackVersion > 0) { /not AOD mode, let's see if it is PicoTrack/
b16bb001	424	AliPicoTrack ptrack = dynamic_cast<AliPicoTrack>(track);
b16bb001	425	if (ptrack) {
cd6431de	426	if (fPicoTrackVersion >= 3)
	427	type = ptrack->GetTrackType();
	428	else
	429	type = ptrack->GetLabel();
b16bb001	430	isEmc = ptrack->IsEMCAL();
	431	}
	432	}
	433
7030f36f	434	Int_t label = 0;
	435
	436	if (fOutMCParticles) {
	437	label = track->GetLabel();
	438
	439	if (label == 0)
	440	label = 99999;
	441	}
	442
ca5c29fa	443	AliDebug(3, Form("Embedding track with pT = %f, eta = %f, phi = %f, label = %d", track->Pt(), track->Eta(), track->Phi(), label));
7030f36f	444	AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc, label);
b16bb001	445	}
	446	}
	447	}
	448
787a3c4f	449	if (fOutClusters) {
b16bb001	450
787a3c4f	451	if (fCopyArray && fClusters)
b16bb001	452	CopyClusters();
	453
	454	if (fAODClusters) {
	455	for (Int_t i = 0; i < fAODClusters->GetEntriesFast(); i++) {
	456	AliVCluster clus = static_cast<AliVCluster>(fAODClusters->At(i));
	457	if (!clus) {
	458	AliError(Form("Could not find cluster %d in branch %s of tree %s!", i, fAODClusName.Data(), fAODTreeName.Data()));
	459	continue;
	460	}
	461	TLorentzVector vect;
	462	Double_t vert[3] = {0,0,0};
	463	clus->GetMomentum(vect,vert);
787a3c4f	464	AddCluster(clus->E(), vect.Eta(), vect.Phi(), clus->GetLabel());
b16bb001	465	}
	466	}
	467	}
	468
787a3c4f	469	if (fOutCaloCells) {
	470
	471	Int_t totalCells = 0;
b16bb001	472
787a3c4f	473	if (fCaloCells)
787a3c4f	474	totalCells += fCaloCells->GetNumberOfCells();
b16bb001	475	if (fAODCaloCells)
	476	totalCells += fAODCaloCells->GetNumberOfCells();
	477
	478	SetNumberOfOutCells(totalCells);
	479
787a3c4f	480	if (fCaloCells)
787a3c4f	481	CopyCells();
b16bb001	482
	483	if (fAODCaloCells) {
	484	for (Short_t i = 0; i < fAODCaloCells->GetNumberOfCells(); i++) {
60d77596	485	Int_t mclabel = 0;
b16bb001	486	Double_t efrac = 0.;
	487	Double_t time = -1;
	488	Short_t cellNum = -1;
	489	Double_t amp = -1;
	490
	491	fAODCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
	492	AliDebug(3,Form("Adding cell with amplitude %f, absolute ID %d, time %f", amp, cellNum, time));
787a3c4f	493	AddCell(amp, cellNum, time, mclabel);
b16bb001	494	}
b16bb001	495	}
cd6431de	496
cd6431de	497	AliDebug(2,Form("Added cells = %d, total cells = %d", fAddedCells, totalCells));
b16bb001	498	}
b16bb001	499	}