[u/mrichter/AliRoot.git] / PWG2 / EBYE / AliAnalysisTaskBF.cxx

#include "TChain.h"
#include "TList.h"
#include "TCanvas.h"
#include "TLorentzVector.h"
#include "TGraphErrors.h"
#include "TH1F.h"

#include "AliAnalysisTaskSE.h"
#include "AliAnalysisManager.h"

#include "AliESDVertex.h"
#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
#include "AliAODEvent.h"
#include "AliAODTrack.h"
#include "AliAODInputHandler.h"
#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
#include "AliStack.h"
#include "AliESDtrackCuts.h"

#include "AliBalance.h"

#include "AliAnalysisTaskBF.h"

// Analysis task for the BF code
// Authors: Panos Cristakoglou@cern.ch

ClassImp(AliAnalysisTaskBF)

//________________________________________________________________________
AliAnalysisTaskBF::AliAnalysisTaskBF(const char *name) 
: AliAnalysisTaskSE(name), 
  fBalance(0),
  fList(0),
  fHistEventStats(0),
  fHistVx(0),
  fHistVy(0),
  fHistVz(0),
  fESDtrackCuts(0),
  fUseOfflineTrigger(kFALSE),
  fVxMax(0.3),
  fVyMax(0.3),
  fVzMax(10.) {
  // Constructor

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 writes into a TH1 container
  DefineOutput(1, AliBalance::Class());
  DefineOutput(2, TList::Class());
}

//________________________________________________________________________
void AliAnalysisTaskBF::UserCreateOutputObjects() {
  // Create histograms
  // Called once
  if(!fBalance) {
    fBalance = new AliBalance();
    fBalance->SetAnalysisLevel("ESD");
    fBalance->SetAnalysisType(1);
    fBalance->SetNumberOfBins(18);
    fBalance->SetInterval(-0.9,0.9);
  }

  //QA list
  fList = new TList();
  fList->SetName("listQA");
  fList->SetOwner();

  //Event stats.
  TString gCutName[4] = {"Total","Offline trigger",
                         "Vertex","Analyzed"};
  fHistEventStats = new TH1F("fHistEventStats",
                             "Event statistics;;N_{events}",
                             4,0.5,4.5);
  for(Int_t i = 1; i <= 4; i++)
    fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
  fList->Add(fHistEventStats);

  //Vertex distributions
  fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);
  fList->Add(fHistVx);
  fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);
  fList->Add(fHistVy);
  fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);
  fList->Add(fHistVz);

  if(fESDtrackCuts) fList->Add(fESDtrackCuts);

  // Post output data.
  PostData(1, fBalance);
  PostData(2, fList);
}

//________________________________________________________________________
void AliAnalysisTaskBF::UserExec(Option_t *) {
  // Main loop
  // Called for each event
  TString gAnalysisLevel = fBalance->GetAnalysisLevel();

  TObjArray *array = new TObjArray();
  
  //ESD analysis
  if(gAnalysisLevel == "ESD") {
    AliESDEvent* gESD = dynamic_cast<AliESDEvent*>(InputEvent()); // from TaskSE
    if (!gESD) {
      Printf("ERROR: gESD not available");
      return;
    }

    fHistEventStats->Fill(1); //all events
    Bool_t isSelected = kTRUE;
    if(fUseOfflineTrigger)
      isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
    if(isSelected) {
      fHistEventStats->Fill(2); //triggered events

      const AliESDVertex *vertex = gESD->GetPrimaryVertex();
      if(vertex) {
	if(vertex->GetNContributors() > 0) {
	  if(vertex->GetZRes() != 0) {
	    fHistEventStats->Fill(3); //events with a proper vertex
	    if(TMath::Abs(vertex->GetXv()) < fVxMax) {
	      if(TMath::Abs(vertex->GetYv()) < fVyMax) {
		if(TMath::Abs(vertex->GetZv()) < fVzMax) {
		  fHistEventStats->Fill(4); //analayzed events
		  fHistVx->Fill(vertex->GetXv());
		  fHistVy->Fill(vertex->GetYv());
		  fHistVz->Fill(vertex->GetZv());

		  //Printf("There are %d tracks in this event", gESD->GetNumberOfTracks());
		  for (Int_t iTracks = 0; iTracks < gESD->GetNumberOfTracks(); iTracks++) {
		    AliESDtrack* track = gESD->GetTrack(iTracks);
		    if (!track) {
		      Printf("ERROR: Could not receive track %d", iTracks);
		      continue;
		    }
		    
		    //ESD track cuts
		    if(fESDtrackCuts) 
		      if(!fESDtrackCuts->AcceptTrack(track)) continue;
		    array->Add(track);
		  } //track loop
		}//Vz cut
	      }//Vy cut
	    }//Vx cut
	  }//proper vertex resolution
	}//proper number of contributors
      }//vertex object valid
    }//triggered event 
  }//ESD analysis
  //AOD analysis
  else if(gAnalysisLevel == "AOD") {
    AliAODEvent* gAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // from TaskSE
    if(!gAOD) {
      Printf("ERROR: gAOD not available");
      return;
    }

    Printf("There are %d tracks in this event", gAOD->GetNumberOfTracks());
    for (Int_t iTracks = 0; iTracks < gAOD->GetNumberOfTracks(); iTracks++) {
      AliAODTrack* track = gAOD->GetTrack(iTracks);
      if (!track) {
	Printf("ERROR: Could not receive track %d", iTracks);
	continue;
      }
      array->Add(track);
    } //track loop
  }//AOD analysis
  //MC analysis
  else if(gAnalysisLevel == "MC") {

    AliMCEvent*  mcEvent = MCEvent(); 
    if (!mcEvent) {
      Printf("ERROR: mcEvent not available");
      return;
    }
    
    Printf("There are %d tracks in this event", mcEvent->GetNumberOfPrimaries());
    for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfPrimaries(); iTracks++) {
      AliMCParticle* track = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(iTracks));
      if (!track) {
	Printf("ERROR: Could not receive particle %d", iTracks);
	continue;
      }
      array->Add(track);
    } //track loop
  }//MC analysis

  fBalance->CalculateBalance(array);
  
  delete array;
  
}      

//________________________________________________________________________
void AliAnalysisTaskBF::Terminate(Option_t *) {
  // Draw result to the screen
  // Called once at the end of the query
  fBalance = dynamic_cast<AliBalance*> (GetOutputData(1));
  if (!fBalance) {
    Printf("ERROR: fBalance not available");
    return;
  }
  
  TGraphErrors *gr = fBalance->DrawBalance();
  gr->SetMarkerStyle(20);
  gr->Draw("AP");

  fBalance->PrintResults();
}
Commit	Line	Data
5c33329d	1	#include "TChain.h"
2cd42194	2	#include "TList.h"
5c33329d	3	#include "TCanvas.h"
	4	#include "TLorentzVector.h"
	5	#include "TGraphErrors.h"
2cd42194	6	#include "TH1F.h"
5c33329d	7
2cd42194	8	#include "AliAnalysisTaskSE.h"
5c33329d	9	#include "AliAnalysisManager.h"
5c33329d	10
2cd42194	11	#include "AliESDVertex.h"
5c33329d	12	#include "AliESDEvent.h"
5c33329d	13	#include "AliESDInputHandler.h"
cd54a838	14	#include "AliAODEvent.h"
6c60d885	15	#include "AliAODTrack.h"
cd54a838	16	#include "AliAODInputHandler.h"
	17	#include "AliMCEventHandler.h"
	18	#include "AliMCEvent.h"
	19	#include "AliStack.h"
5115b492	20	#include "AliESDtrackCuts.h"
5c33329d	21
	22	#include "AliBalance.h"
	23
	24	#include "AliAnalysisTaskBF.h"
	25
	26	// Analysis task for the BF code
	27	// Authors: Panos Cristakoglou@cern.ch
	28
	29	ClassImp(AliAnalysisTaskBF)
	30
	31	//________________________________________________________________________
	32	AliAnalysisTaskBF::AliAnalysisTaskBF(const char *name)
09bb7bf4	33	: AliAnalysisTaskSE(name),
	34	fBalance(0),
	35	fList(0),
	36	fHistEventStats(0),
	37	fHistVx(0),
	38	fHistVy(0),
	39	fHistVz(0),
	40	fESDtrackCuts(0),
	41	fUseOfflineTrigger(kFALSE),
	42	fVxMax(0.3),
	43	fVyMax(0.3),
	44	fVzMax(10.) {
5c33329d	45	// Constructor
	46
	47	// Define input and output slots here
	48	// Input slot #0 works with a TChain
	49	DefineInput(0, TChain::Class());
	50	// Output slot #0 writes into a TH1 container
2cd42194	51	DefineOutput(1, AliBalance::Class());
2cd42194	52	DefineOutput(2, TList::Class());
5c33329d	53	}
	54
	55	//________________________________________________________________________
2cd42194	56	void AliAnalysisTaskBF::UserCreateOutputObjects() {
2cd42194	57	// Create histograms
5c33329d	58	// Called once
2cd42194	59	if(!fBalance) {
	60	fBalance = new AliBalance();
	61	fBalance->SetAnalysisLevel("ESD");
	62	fBalance->SetAnalysisType(1);
	63	fBalance->SetNumberOfBins(18);
	64	fBalance->SetInterval(-0.9,0.9);
5c33329d	65	}
5c33329d	66
09bb7bf4	67	//QA list
2cd42194	68	fList = new TList();
2cd42194	69	fList->SetName("listQA");
0b01238f	70	fList->SetOwner();
2cd42194	71
09bb7bf4	72	//Event stats.
2cd42194	73	TString gCutName[4] = {"Total","Offline trigger",
	74	"Vertex","Analyzed"};
	75	fHistEventStats = new TH1F("fHistEventStats",
	76	"Event statistics;;N_{events}",
	77	4,0.5,4.5);
	78	for(Int_t i = 1; i <= 4; i++)
	79	fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
	80	fList->Add(fHistEventStats);
	81
09bb7bf4	82	//Vertex distributions
	83	fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);
	84	fList->Add(fHistVx);
	85	fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);
	86	fList->Add(fHistVy);
	87	fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);
	88	fList->Add(fHistVz);
	89
5115b492	90	if(fESDtrackCuts) fList->Add(fESDtrackCuts);
5115b492	91
2cd42194	92	// Post output data.
	93	PostData(1, fBalance);
	94	PostData(2, fList);
5c33329d	95	}
	96
	97	//________________________________________________________________________
2cd42194	98	void AliAnalysisTaskBF::UserExec(Option_t *) {
5c33329d	99	// Main loop
5c33329d	100	// Called for each event
cd54a838	101	TString gAnalysisLevel = fBalance->GetAnalysisLevel();
5c33329d	102
5c33329d	103	TObjArray *array = new TObjArray();
cd54a838	104
	105	//ESD analysis
	106	if(gAnalysisLevel == "ESD") {
2cd42194	107	AliESDEvent* gESD = dynamic_cast<AliESDEvent*>(InputEvent()); // from TaskSE
	108	if (!gESD) {
	109	Printf("ERROR: gESD not available");
cd54a838	110	return;
5c33329d	111	}
2cd42194	112
2cd42194	113	fHistEventStats->Fill(1); //all events
09bb7bf4	114	Bool_t isSelected = kTRUE;
	115	if(fUseOfflineTrigger)
	116	isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
2cd42194	117	if(isSelected) {
	118	fHistEventStats->Fill(2); //triggered events
	119
	120	const AliESDVertex *vertex = gESD->GetPrimaryVertex();
	121	if(vertex) {
7f94638c	122	if(vertex->GetNContributors() > 0) {
	123	if(vertex->GetZRes() != 0) {
	124	fHistEventStats->Fill(3); //events with a proper vertex
	125	if(TMath::Abs(vertex->GetXv()) < fVxMax) {
	126	if(TMath::Abs(vertex->GetYv()) < fVyMax) {
	127	if(TMath::Abs(vertex->GetZv()) < fVzMax) {
	128	fHistEventStats->Fill(4); //analayzed events
09bb7bf4	129	fHistVx->Fill(vertex->GetXv());
	130	fHistVy->Fill(vertex->GetYv());
	131	fHistVz->Fill(vertex->GetZv());
	132
	133	//Printf("There are %d tracks in this event", gESD->GetNumberOfTracks());
7f94638c	134	for (Int_t iTracks = 0; iTracks < gESD->GetNumberOfTracks(); iTracks++) {
	135	AliESDtrack* track = gESD->GetTrack(iTracks);
	136	if (!track) {
	137	Printf("ERROR: Could not receive track %d", iTracks);
	138	continue;
	139	}
	140
	141	//ESD track cuts
	142	if(fESDtrackCuts)
	143	if(!fESDtrackCuts->AcceptTrack(track)) continue;
	144	array->Add(track);
	145	} //track loop
	146	}//Vz cut
	147	}//Vy cut
	148	}//Vx cut
	149	}//proper vertex resolution
	150	}//proper number of contributors
2cd42194	151	}//vertex object valid
2cd42194	152	}//triggered event
cd54a838	153	}//ESD analysis
	154	//AOD analysis
	155	else if(gAnalysisLevel == "AOD") {
2cd42194	156	AliAODEvent* gAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // from TaskSE
	157	if(!gAOD) {
	158	Printf("ERROR: gAOD not available");
cd54a838	159	return;
cd54a838	160	}
2cd42194	161
	162	Printf("There are %d tracks in this event", gAOD->GetNumberOfTracks());
	163	for (Int_t iTracks = 0; iTracks < gAOD->GetNumberOfTracks(); iTracks++) {
	164	AliAODTrack* track = gAOD->GetTrack(iTracks);
cd54a838	165	if (!track) {
	166	Printf("ERROR: Could not receive track %d", iTracks);
	167	continue;
	168	}
	169	array->Add(track);
	170	} //track loop
	171	}//AOD analysis
	172	//MC analysis
	173	else if(gAnalysisLevel == "MC") {
2cd42194	174
	175	AliMCEvent* mcEvent = MCEvent();
	176	if (!mcEvent) {
	177	Printf("ERROR: mcEvent not available");
cd54a838	178	return;
	179	}
	180
2cd42194	181	Printf("There are %d tracks in this event", mcEvent->GetNumberOfPrimaries());
	182	for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfPrimaries(); iTracks++) {
	183	AliMCParticle* track = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(iTracks));
cd54a838	184	if (!track) {
	185	Printf("ERROR: Could not receive particle %d", iTracks);
	186	continue;
	187	}
	188	array->Add(track);
	189	} //track loop
	190	}//MC analysis
	191
5c33329d	192	fBalance->CalculateBalance(array);
	193
	194	delete array;
2cd42194	195
5c33329d	196	}
	197
	198	//________________________________________________________________________
	199	void AliAnalysisTaskBF::Terminate(Option_t *) {
	200	// Draw result to the screen
	201	// Called once at the end of the query
2cd42194	202	fBalance = dynamic_cast<AliBalance*> (GetOutputData(1));
5c33329d	203	if (!fBalance) {
	204	Printf("ERROR: fBalance not available");
	205	return;
	206	}
	207
	208	TGraphErrors *gr = fBalance->DrawBalance();
	209	gr->SetMarkerStyle(20);
	210	gr->Draw("AP");
	211
	212	fBalance->PrintResults();
	213	}