[u/mrichter/AliRoot.git] / PWGPP / VZERO / AliAnaVZEROTrigger.cxx

#include <stdio.h>
#include <stdlib.h>

#include "TH1F.h"
#include "TH2F.h"
#include "TChain.h"
#include "TFile.h"
#include "TParameter.h"

#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"

#include "AliESDEvent.h"
#include "AliESDfriend.h"
#include "AliESDInputHandler.h"
#include "AliAnaVZEROTrigger.h"
#include "AliCentrality.h"

// VZERO includes
#include "AliESDVZERO.h"

ClassImp(AliAnaVZEROTrigger)

AliAnaVZEROTrigger::AliAnaVZEROTrigger() 
  : AliAnalysisTaskSE("AliAnaVZEROTrigger"), fESD(0), fOutputList(0),
  fMinThr(400.),
  fMaxThr(14000.),
  fRatio(2.0),
  fNThr(30),
  fMBTrigName("CPBI"),
  fUsePhysSel(kFALSE),
  fV0Percent(0),
  fV0PercentAll(0),
  fZvtx(0),
  fXYvtx(0),
  fV0Mult1d(0),
  fV0Charge2d(0),
  fV0Charge2dPercent(0),
  fV0Charge2dAll(0),
  fV0PercentBins(0),
  fV0PercentBinsAll(0),
  fV0Cent(0),
  fV0CentAll(0),
  fV0SemiCent(0),
  fV0SemiCentAll(0),
  fV0CentHw(0),
  fV0CentHwAll(0),
  fV0SemiCentHw(0),
  fV0SemiCentHwAll(0),
  fV0CentTr(0),
  fV0CentTrAll(0),
  fV0SemiCentTr(0),
  fV0SemiCentTrAll(0),
  fV0Percent63(0),
  fV0Percent63All(0),
  fV0MultAll(0),
  fV0Mult63(0),
  fV0CentVsMult(0),
  fV0CentVsCharge(0),
  fV0CentVsTrCharge(0),
  fV0TrChargeVsChargeA(0),
  fV0TrChargeVsChargeC(0),
  fAdcPmt(0),
  fAdcPmt0(0),
  fAdcPmt1(0),
  fAdcPmt2(0),
  fAdcPmt3(0)
{
  // Constructor
  // Init default thr
  fCentCuts[0] = 4152;
  fCentCuts[1] = 8637;
  fSemiCentCuts[0] = 617;
  fSemiCentCuts[1] = 1283;

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 id reserved by the base class for AOD
  // Output slot #1 writes into a TH1 container
  DefineOutput(1, TList::Class());
}

//________________________________________________________________________
AliAnaVZEROTrigger::AliAnaVZEROTrigger(const char *name) 
  : AliAnalysisTaskSE(name), fESD(0), fOutputList(0),
  fMinThr(400.),
  fMaxThr(14000.),
  fRatio(2.0),
  fNThr(30),
  fMBTrigName("CPBI"),
  fUsePhysSel(kFALSE),
  fV0Percent(0),
  fV0PercentAll(0),
  fZvtx(0),
  fXYvtx(0),
  fV0Mult1d(0),
  fV0Charge2d(0),
  fV0Charge2dPercent(0),
  fV0Charge2dAll(0),
  fV0PercentBins(0),
  fV0PercentBinsAll(0),
  fV0Cent(0),
  fV0CentAll(0),
  fV0SemiCent(0),
  fV0SemiCentAll(0),
  fV0CentHw(0),
  fV0CentHwAll(0),
  fV0SemiCentHw(0),
  fV0SemiCentHwAll(0),
  fV0CentTr(0),
  fV0CentTrAll(0),
  fV0SemiCentTr(0),
  fV0SemiCentTrAll(0),
  fV0Percent63(0),
  fV0Percent63All(0),
  fV0MultAll(0),
  fV0Mult63(0),
  fV0CentVsMult(0),
  fV0CentVsCharge(0),
  fV0CentVsTrCharge(0),
  fV0TrChargeVsChargeA(0),
  fV0TrChargeVsChargeC(0),
  fAdcPmt(0),
  fAdcPmt0(0),
  fAdcPmt1(0),
  fAdcPmt2(0),
  fAdcPmt3(0)
{
  // Constructor
  // Init default thr
  fCentCuts[0] = 4152;
  fCentCuts[1] = 8637;
  fSemiCentCuts[0] = 617;
  fSemiCentCuts[1] = 1283;

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 id reserved by the base class for AOD
  // Output slot #1 writes into a TH1 container
  DefineOutput(1, TList::Class());
}

//________________________________________________________________________
void AliAnaVZEROTrigger::Setup(const char *filename)
{
  // Open the config file and
  // read the parameters
  FILE *fin = fopen(filename,"r");
  if (!fin) {
    AliInfo(Form("File %s is not found running with the default parameters.",filename));
  }
  else {
    Int_t res = fscanf(fin,"%f %f %f %d %f %f %f %f",
		       &fMinThr,&fMaxThr,&fRatio,&fNThr,
		       &fSemiCentCuts[0],&fSemiCentCuts[1],&fCentCuts[0],&fCentCuts[1]);
    if(res!=8) {
      AliFatal("Failed to get values from the config file.\n");
    }
    fclose(fin);
  }

  AliInfo(Form("MinThr=%.1f MaxThr=%.1f Ratio=%.4f NThr=%d",fMinThr,fMaxThr,fRatio,fNThr));
  AliInfo(Form("CTA1=%.1f CTC1=%.1f CTA2=%.1f CTC2=%.1f",fSemiCentCuts[0],fSemiCentCuts[1],fCentCuts[0],fCentCuts[1]));
}

//________________________________________________________________________
void AliAnaVZEROTrigger::UserCreateOutputObjects()
{
  // Create histograms
  // Called once

  fOutputList = new TList();
  fOutputList->SetOwner(kTRUE);

  fV0Percent    = new TH1F("fV0Percent","Centrality percentile based on v0 mult",400,0,100);
  fOutputList->Add(fV0Percent);
  fV0PercentAll    = new TH1F("fV0PercentAll","Centrality percentile based on v0 mult (no event selection)",400,0,100);
  fOutputList->Add(fV0PercentAll);
  fZvtx    = new TH1F("fZvtx","",500,-20,20);
  fOutputList->Add(fZvtx);
  fXYvtx    = new TH2F("fXYvtx","",250,-0.5,0.5,250,-0.5,0.5);
  fOutputList->Add(fXYvtx);
  fV0Mult1d     = new TH1F("fV0Mult1d","Total v0 mult",500,0,25000);
  fOutputList->Add(fV0Mult1d);

  fV0Charge2d     = new TH2F("fV0Charge2d","V0C vs V0A charge",125,0,25000,125,0,25000);
  fOutputList->Add(fV0Charge2d);
  fV0Charge2dPercent     = new TH2F("fV0Charge2dPercent","V0C vs V0A charge vs centrality percentile",125,0,25000,125,0,25000);
  fOutputList->Add(fV0Charge2dPercent);
  fV0Charge2dAll  = new TH2F("fV0Charge2dAll","V0C vs V0A charge (all events)",125,0,25000,125,0,25000);
  fOutputList->Add(fV0Charge2dAll);

  fV0PercentBins = new TH1F*[fNThr];
  fV0PercentBinsAll = new TH1F*[fNThr];
  for(Int_t j = 0; j < fNThr; ++j) {
    fV0PercentBins[j] = new TH1F(Form("fV0PercentBins_%d",j),"Centrality percentile with V0 charge thresholds",400,0,100);
    fOutputList->Add(fV0PercentBins[j]);
    fV0PercentBinsAll[j] = new TH1F(Form("fV0PercentBinsAll_%d",j),"Centrality percentile with V0 charge thresholds (no event selection)",100,0,100);
    fOutputList->Add(fV0PercentBinsAll[j]);
  }

  fV0Cent = new TH1F("fV0Cent","Centrality percentile with custom V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0Cent);
  fV0CentAll = new TH1F("fV0CentAll","Centrality percentile with custom V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0CentAll);
  fV0SemiCent = new TH1F("fV0SemiCent","Centrality percentile with custom V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0SemiCent);
  fV0SemiCentAll = new TH1F("fV0SemiCentAll","Centrality percentile with custom V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0SemiCentAll);

  fV0CentHw = new TH1F("fV0CentHw","Centrality percentile with hardware V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0CentHw);
  fV0CentHwAll = new TH1F("fV0CentHwAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0CentHwAll);
  fV0SemiCentHw = new TH1F("fV0SemiCentHw","Centrality percentile with hardware V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0SemiCentHw);
  fV0SemiCentHwAll = new TH1F("fV0SemiCentHwAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0SemiCentHwAll);

  fV0CentTr = new TH1F("fV0CentTr","Centrality percentile with hardware V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0CentTr);
  fV0CentTrAll = new TH1F("fV0CentTrAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0CentTrAll);
  fV0SemiCentTr = new TH1F("fV0SemiCentTr","Centrality percentile with hardware V0 charge thresholds",400,0,100);
  fOutputList->Add(fV0SemiCentTr);
  fV0SemiCentTrAll = new TH1F("fV0SemiCentTrAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
  fOutputList->Add(fV0SemiCentTrAll);

  fV0Percent63    = new TH1F("fV0Percent63","Centrality percentile based on v0 mult",400,0,100);
  fOutputList->Add(fV0Percent63);
  fV0Percent63All    = new TH1F("fV0Percent63All","Centrality percentile based on v0 mult (no event selection)",400,0,100);
  fOutputList->Add(fV0Percent63All);

  fV0MultAll = new TH2F("fV0MultAll","",250,0,25000,250,0,25000);
  fOutputList->Add(fV0MultAll);
  fV0Mult63 = new TH2F("fV0Mult63","",250,0,25000,250,0,25000);
  fOutputList->Add(fV0Mult63);

  fV0CentVsMult = new TH2F("fV0CentVsMult","",250,0,25000,200,0,100);
  fOutputList->Add(fV0CentVsMult);
  fV0CentVsCharge = new TH2F("fV0CentVsCharge","",250,0,100000,200,0,100);
  fOutputList->Add(fV0CentVsCharge); 
  fV0CentVsTrCharge = new TH2F("fV0CentVsTrCharge","",250,0,50000,200,0,100);
  fOutputList->Add(fV0CentVsTrCharge);
  fV0TrChargeVsChargeA = new TH2F("fV0TrChargeVsChargeA","",250,0,75000,250,0,37500);
  fOutputList->Add(fV0TrChargeVsChargeA);
  fV0TrChargeVsChargeC = new TH2F("fV0TrChargeVsChargeC","",250,0,75000,250,0,37500);
  fOutputList->Add(fV0TrChargeVsChargeC);

  fAdcPmt = new TH2F("fAdcPmt","",64,-0.5,63.5,400,0,4000);
  fOutputList->Add(fAdcPmt);
  fAdcPmt0 = new TH2F("fAdcPmt0","",64,-0.5,63.5,400,0,4000);
  fOutputList->Add(fAdcPmt0);
  fAdcPmt1 = new TH2F("fAdcPmt1","",64,-0.5,63.5,400,0,4000);
  fOutputList->Add(fAdcPmt1);
  fAdcPmt2 = new TH2F("fAdcPmt2","",64,-0.5,63.5,400,0,4000);
  fOutputList->Add(fAdcPmt2);
  fAdcPmt3 = new TH2F("fAdcPmt3","",64,-0.5,63.5,400,0,4000);
  fOutputList->Add(fAdcPmt3);

  PostData(1, fOutputList);
 }
//________________________________________________________________________
void AliAnaVZEROTrigger::Init()
{
  // Nothing here
  // ....
}

//________________________________________________________________________
void AliAnaVZEROTrigger::UserExec(Option_t *) 
{
  // Main loop
  // Called for each event


  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
  if (!fESD) {
    printf("ERROR: fESD not available\n");
    return;
  }

  AliESDVZERO* esdV0 = fESD->GetVZEROData();
  if (!esdV0) {
    Printf("ERROR: esd V0  not available");
    return;
  }

  // Trigger
  TString trigStr(fESD->GetFiredTriggerClasses());
  if (!trigStr.Contains("-B-")) return;
  if (!trigStr.Contains(fMBTrigName.Data())) return;

  // Count V0 flags
  Int_t nV0A = 0;
  Int_t nV0C = 0;
  for(Int_t i = 0; i < 32; ++i) {
    if (esdV0->GetBBFlag(i)) nV0C++;
    if (esdV0->GetBBFlag(i+32)) nV0A++;
  }

  // Phys sel
  Bool_t goodEvent = kTRUE;
  Bool_t isSelected;
  if (fUsePhysSel)
    isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
  else
    isSelected = ((esdV0->GetV0ADecision()==1) && (esdV0->GetV0CDecision()==1));

  if (!isSelected) goodEvent = kFALSE;

  const AliESDVertex *primaryVtx = fESD->GetPrimaryVertexSPD();
  if (!primaryVtx) goodEvent = kFALSE;
  if (!primaryVtx->GetStatus()) goodEvent = kFALSE;
  Double_t tPrimaryVtxPosition[3];
  primaryVtx->GetXYZ(tPrimaryVtxPosition);
  if (goodEvent) {
    fZvtx->Fill(tPrimaryVtxPosition[2]);
    fXYvtx->Fill(tPrimaryVtxPosition[0],tPrimaryVtxPosition[1]);
  }
  if (TMath::Abs(tPrimaryVtxPosition[2]) > 10.0) goodEvent = kFALSE;

  if (goodEvent) {
    for(Int_t i = 0; i < 64; ++i) {
      fAdcPmt->Fill(i,esdV0->GetAdc(i));
      if (tPrimaryVtxPosition[2]>5.0) fAdcPmt0->Fill(i,esdV0->GetAdc(i));
      if (tPrimaryVtxPosition[2]>0.0 && tPrimaryVtxPosition[2]<5.0) fAdcPmt1->Fill(i,esdV0->GetAdc(i));
      if (tPrimaryVtxPosition[2]>-5.0 && tPrimaryVtxPosition[2]<0.0) fAdcPmt2->Fill(i,esdV0->GetAdc(i));
      if (tPrimaryVtxPosition[2]<-5.0) fAdcPmt3->Fill(i,esdV0->GetAdc(i));
    }
  }

  UShort_t chargeA = esdV0->GetTriggerChargeA();
  UShort_t chargeC = esdV0->GetTriggerChargeC();

  Float_t offChargeA = 0;
  Float_t offChargeC = 0;
  for(Int_t i = 0; i < 64; ++i) {
    if (i < 32) offChargeC += esdV0->GetAdc(i);
    else offChargeA += esdV0->GetAdc(i);
  }

  if (goodEvent) {
    fV0TrChargeVsChargeA->Fill(offChargeA,chargeA);
    fV0TrChargeVsChargeC->Fill(offChargeC,chargeC);
  }

  Float_t percentile = 0;
  AliCentrality *centrality = fESD->GetCentrality();
  //  if (centrality->GetQuality()) goodEvent = kFALSE;
  percentile = centrality->GetCentralityPercentile("V0M");
  //  if (percentile < 0) percentile = 0;

  if (goodEvent) {
    fV0CentVsMult->Fill(esdV0->GetMTotV0A()+esdV0->GetMTotV0C(),percentile);
    fV0CentVsCharge->Fill(offChargeA+offChargeC,percentile);
    fV0CentVsTrCharge->Fill(chargeA+chargeC,percentile);
  }

  fV0MultAll->Fill(esdV0->GetMTotV0A(),esdV0->GetMTotV0C());
  if ((nV0A+nV0C)>=63)
    fV0Mult63->Fill(esdV0->GetMTotV0A(),esdV0->GetMTotV0C());

  fV0PercentAll->Fill(percentile);
  if (goodEvent) fV0Percent->Fill(percentile);

  if ((nV0A+nV0C)>=63) {
    fV0Percent63All->Fill(percentile);
    if (goodEvent) fV0Percent63->Fill(percentile);
  }

  Float_t multV0 = esdV0->GetMTotV0A()+esdV0->GetMTotV0C();
  if (goodEvent) fV0Mult1d->Fill(multV0);

  fV0Charge2dAll->Fill(chargeA,chargeC);
  if (goodEvent) {
    fV0Charge2d->Fill(chargeA,chargeC);
    fV0Charge2dPercent->Fill(chargeA,chargeC,percentile);
  }

  for(Int_t j = 0; j < fNThr; ++j) {
    Float_t thrA = GetThrA(j);
    Float_t thrC = GetThrC(j);
    if (chargeA >= thrA && chargeC >= thrC) {
      fV0PercentBinsAll[j]->Fill(percentile);
      if (goodEvent) fV0PercentBins[j]->Fill(percentile);
    }
  }

  if (chargeA >= fCentCuts[0] && chargeC >= fCentCuts[1]) {
    fV0CentAll->Fill(percentile);
    if (goodEvent) fV0Cent->Fill(percentile);
  }
  if (chargeA >= fSemiCentCuts[0] && chargeC >= fSemiCentCuts[1]) {
    fV0SemiCentAll->Fill(percentile);
    if (goodEvent) fV0SemiCent->Fill(percentile);
  }

  if (esdV0->GetTriggerBits() & (1<<8)) {
    fV0CentHwAll->Fill(percentile);
    if (goodEvent) fV0CentHw->Fill(percentile);
  }
  if (esdV0->GetTriggerBits() & (1<<6)) {
    fV0SemiCentHwAll->Fill(percentile);
    if (goodEvent) fV0SemiCentHw->Fill(percentile);
  }

  if (trigStr.Contains("CVHN")) {
    fV0CentTrAll->Fill(percentile);
    if (goodEvent) fV0CentTr->Fill(percentile);
  }
  if (trigStr.Contains("CVLN")) {
    fV0SemiCentTrAll->Fill(percentile);
    if (goodEvent) fV0SemiCentTr->Fill(percentile);
  }

  PostData(1, fOutputList);
}      

//________________________________________________________________________
void AliAnaVZEROTrigger::Terminate(Option_t *) 
{
  // Check the output list and store output config file
  // Called once at the end of the query

  fOutputList = dynamic_cast<TList*> (GetOutputData(1));
  if (!fOutputList) {
    printf("ERROR: Output list not available\n");
    return;
  }

  printf("\n\n");
  FILE *fout=fopen("./trigger.txt","w");
  if (!fout) {
    printf("Failed to open local result file\n");
    return;
  }
  printf("%.0f %.0f %f %d %.0f %.0f %.0f %.0f\n\n",
	 GetMinThr(),
	 GetMaxThr(),
	 GetRatio(),
	 GetNThr(),
	 GetSemiCentCuts(0),GetSemiCentCuts(1),
	 GetCentCuts(0),GetCentCuts(1));

  fprintf(fout,"%.0f %.0f %f %d %.0f %.0f %.0f %.0f\n\n",
	  GetMinThr(),
	  GetMaxThr(),
	  GetRatio(),
	  GetNThr(),
	  GetSemiCentCuts(0),GetSemiCentCuts(1),
	  GetCentCuts(0),GetCentCuts(1));
  fclose(fout);
}

//________________________________________________________________________
Float_t AliAnaVZEROTrigger::GetThrA(Int_t j) const
{
  // Get the threshold on A side with
  // index i
  if (j<0 || j>= fNThr) return 0;

  return (fMinThr + ((Float_t)j)*(fMaxThr-fMinThr)/((Float_t)fNThr-1.));
}

//________________________________________________________________________
Float_t AliAnaVZEROTrigger::GetThrC(Int_t j) const
{
  // Get the threshold on C side with
  // index i
  if (j<0 || j>= fNThr) return 0;

  Float_t thrA = (fMinThr + ((Float_t)j)*(fMaxThr-fMinThr)/((Float_t)fNThr-1.));
  return (thrA*fRatio);
}
Commit	Line	Data
c014f07d	1	#include <stdio.h>
	2	#include <stdlib.h>
	3
	4	#include "TH1F.h"
	5	#include "TH2F.h"
	6	#include "TChain.h"
	7	#include "TFile.h"
	8	#include "TParameter.h"
	9
	10	#include "AliAnalysisTask.h"
	11	#include "AliAnalysisManager.h"
	12
	13	#include "AliESDEvent.h"
	14	#include "AliESDfriend.h"
	15	#include "AliESDInputHandler.h"
	16	#include "AliAnaVZEROTrigger.h"
	17	#include "AliCentrality.h"
	18
	19	// VZERO includes
	20	#include "AliESDVZERO.h"
	21
	22	ClassImp(AliAnaVZEROTrigger)
	23
	24	AliAnaVZEROTrigger::AliAnaVZEROTrigger()
	25	: AliAnalysisTaskSE("AliAnaVZEROTrigger"), fESD(0), fOutputList(0),
e5f60c5a	26	fMinThr(400.),
	27	fMaxThr(14000.),
	28	fRatio(2.0),
	29	fNThr(30),
c014f07d	30	fMBTrigName("CPBI"),
e5f60c5a	31	fUsePhysSel(kFALSE),
c014f07d	32	fV0Percent(0),
	33	fV0PercentAll(0),
	34	fZvtx(0),
e5f60c5a	35	fXYvtx(0),
c014f07d	36	fV0Mult1d(0),
	37	fV0Charge2d(0),
	38	fV0Charge2dPercent(0),
	39	fV0Charge2dAll(0),
	40	fV0PercentBins(0),
	41	fV0PercentBinsAll(0),
	42	fV0Cent(0),
	43	fV0CentAll(0),
	44	fV0SemiCent(0),
	45	fV0SemiCentAll(0),
	46	fV0CentHw(0),
	47	fV0CentHwAll(0),
	48	fV0SemiCentHw(0),
	49	fV0SemiCentHwAll(0),
	50	fV0CentTr(0),
	51	fV0CentTrAll(0),
	52	fV0SemiCentTr(0),
	53	fV0SemiCentTrAll(0),
	54	fV0Percent63(0),
	55	fV0Percent63All(0),
	56	fV0MultAll(0),
e5f60c5a	57	fV0Mult63(0),
	58	fV0CentVsMult(0),
	59	fV0CentVsCharge(0),
	60	fV0CentVsTrCharge(0),
	61	fV0TrChargeVsChargeA(0),
	62	fV0TrChargeVsChargeC(0),
	63	fAdcPmt(0),
	64	fAdcPmt0(0),
	65	fAdcPmt1(0),
	66	fAdcPmt2(0),
	67	fAdcPmt3(0)
c014f07d	68	{
	69	// Constructor
	70	// Init default thr
	71	fCentCuts[0] = 4152;
	72	fCentCuts[1] = 8637;
	73	fSemiCentCuts[0] = 617;
	74	fSemiCentCuts[1] = 1283;
	75
	76	// Define input and output slots here
	77	// Input slot #0 works with a TChain
	78	DefineInput(0, TChain::Class());
	79	// Output slot #0 id reserved by the base class for AOD
	80	// Output slot #1 writes into a TH1 container
	81	DefineOutput(1, TList::Class());
	82	}
	83
	84	//________________________________________________________________________
	85	AliAnaVZEROTrigger::AliAnaVZEROTrigger(const char *name)
	86	: AliAnalysisTaskSE(name), fESD(0), fOutputList(0),
e5f60c5a	87	fMinThr(400.),
	88	fMaxThr(14000.),
	89	fRatio(2.0),
	90	fNThr(30),
c014f07d	91	fMBTrigName("CPBI"),
e5f60c5a	92	fUsePhysSel(kFALSE),
c014f07d	93	fV0Percent(0),
	94	fV0PercentAll(0),
	95	fZvtx(0),
e5f60c5a	96	fXYvtx(0),
c014f07d	97	fV0Mult1d(0),
	98	fV0Charge2d(0),
	99	fV0Charge2dPercent(0),
	100	fV0Charge2dAll(0),
	101	fV0PercentBins(0),
	102	fV0PercentBinsAll(0),
	103	fV0Cent(0),
	104	fV0CentAll(0),
	105	fV0SemiCent(0),
	106	fV0SemiCentAll(0),
	107	fV0CentHw(0),
	108	fV0CentHwAll(0),
	109	fV0SemiCentHw(0),
	110	fV0SemiCentHwAll(0),
	111	fV0CentTr(0),
	112	fV0CentTrAll(0),
	113	fV0SemiCentTr(0),
	114	fV0SemiCentTrAll(0),
	115	fV0Percent63(0),
	116	fV0Percent63All(0),
	117	fV0MultAll(0),
e5f60c5a	118	fV0Mult63(0),
	119	fV0CentVsMult(0),
	120	fV0CentVsCharge(0),
	121	fV0CentVsTrCharge(0),
	122	fV0TrChargeVsChargeA(0),
	123	fV0TrChargeVsChargeC(0),
	124	fAdcPmt(0),
	125	fAdcPmt0(0),
	126	fAdcPmt1(0),
	127	fAdcPmt2(0),
	128	fAdcPmt3(0)
c014f07d	129	{
	130	// Constructor
	131	// Init default thr
	132	fCentCuts[0] = 4152;
	133	fCentCuts[1] = 8637;
	134	fSemiCentCuts[0] = 617;
	135	fSemiCentCuts[1] = 1283;
	136
	137	// Define input and output slots here
	138	// Input slot #0 works with a TChain
	139	DefineInput(0, TChain::Class());
	140	// Output slot #0 id reserved by the base class for AOD
	141	// Output slot #1 writes into a TH1 container
	142	DefineOutput(1, TList::Class());
	143	}
	144
	145	//________________________________________________________________________
	146	void AliAnaVZEROTrigger::Setup(const char *filename)
	147	{
	148	// Open the config file and
	149	// read the parameters
	150	FILE *fin = fopen(filename,"r");
	151	if (!fin) {
	152	AliInfo(Form("File %s is not found running with the default parameters.",filename));
	153	}
	154	else {
	155	Int_t res = fscanf(fin,"%f %f %f %d %f %f %f %f",
	156	&fMinThr,&fMaxThr,&fRatio,&fNThr,
	157	&fSemiCentCuts[0],&fSemiCentCuts[1],&fCentCuts[0],&fCentCuts[1]);
	158	if(res!=8) {
	159	AliFatal("Failed to get values from the config file.\n");
	160	}
	161	fclose(fin);
	162	}
	163
	164	AliInfo(Form("MinThr=%.1f MaxThr=%.1f Ratio=%.4f NThr=%d",fMinThr,fMaxThr,fRatio,fNThr));
	165	AliInfo(Form("CTA1=%.1f CTC1=%.1f CTA2=%.1f CTC2=%.1f",fSemiCentCuts[0],fSemiCentCuts[1],fCentCuts[0],fCentCuts[1]));
	166	}
	167
	168	//________________________________________________________________________
	169	void AliAnaVZEROTrigger::UserCreateOutputObjects()
	170	{
	171	// Create histograms
	172	// Called once
	173
	174	fOutputList = new TList();
	175	fOutputList->SetOwner(kTRUE);
	176
	177	fV0Percent = new TH1F("fV0Percent","Centrality percentile based on v0 mult",400,0,100);
	178	fOutputList->Add(fV0Percent);
	179	fV0PercentAll = new TH1F("fV0PercentAll","Centrality percentile based on v0 mult (no event selection)",400,0,100);
	180	fOutputList->Add(fV0PercentAll);
	181	fZvtx = new TH1F("fZvtx","",500,-20,20);
	182	fOutputList->Add(fZvtx);
e5f60c5a	183	fXYvtx = new TH2F("fXYvtx","",250,-0.5,0.5,250,-0.5,0.5);
e5f60c5a	184	fOutputList->Add(fXYvtx);
c014f07d	185	fV0Mult1d = new TH1F("fV0Mult1d","Total v0 mult",500,0,25000);
	186	fOutputList->Add(fV0Mult1d);
	187
	188	fV0Charge2d = new TH2F("fV0Charge2d","V0C vs V0A charge",125,0,25000,125,0,25000);
	189	fOutputList->Add(fV0Charge2d);
	190	fV0Charge2dPercent = new TH2F("fV0Charge2dPercent","V0C vs V0A charge vs centrality percentile",125,0,25000,125,0,25000);
	191	fOutputList->Add(fV0Charge2dPercent);
	192	fV0Charge2dAll = new TH2F("fV0Charge2dAll","V0C vs V0A charge (all events)",125,0,25000,125,0,25000);
	193	fOutputList->Add(fV0Charge2dAll);
	194
	195	fV0PercentBins = new TH1F*[fNThr];
	196	fV0PercentBinsAll = new TH1F*[fNThr];
	197	for(Int_t j = 0; j < fNThr; ++j) {
	198	fV0PercentBins[j] = new TH1F(Form("fV0PercentBins_%d",j),"Centrality percentile with V0 charge thresholds",400,0,100);
	199	fOutputList->Add(fV0PercentBins[j]);
	200	fV0PercentBinsAll[j] = new TH1F(Form("fV0PercentBinsAll_%d",j),"Centrality percentile with V0 charge thresholds (no event selection)",100,0,100);
	201	fOutputList->Add(fV0PercentBinsAll[j]);
	202	}
	203
	204	fV0Cent = new TH1F("fV0Cent","Centrality percentile with custom V0 charge thresholds",400,0,100);
	205	fOutputList->Add(fV0Cent);
	206	fV0CentAll = new TH1F("fV0CentAll","Centrality percentile with custom V0 charge thresholds (no event selection)",100,0,100);
	207	fOutputList->Add(fV0CentAll);
	208	fV0SemiCent = new TH1F("fV0SemiCent","Centrality percentile with custom V0 charge thresholds",400,0,100);
	209	fOutputList->Add(fV0SemiCent);
	210	fV0SemiCentAll = new TH1F("fV0SemiCentAll","Centrality percentile with custom V0 charge thresholds (no event selection)",100,0,100);
	211	fOutputList->Add(fV0SemiCentAll);
	212
	213	fV0CentHw = new TH1F("fV0CentHw","Centrality percentile with hardware V0 charge thresholds",400,0,100);
	214	fOutputList->Add(fV0CentHw);
	215	fV0CentHwAll = new TH1F("fV0CentHwAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
	216	fOutputList->Add(fV0CentHwAll);
	217	fV0SemiCentHw = new TH1F("fV0SemiCentHw","Centrality percentile with hardware V0 charge thresholds",400,0,100);
	218	fOutputList->Add(fV0SemiCentHw);
	219	fV0SemiCentHwAll = new TH1F("fV0SemiCentHwAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
	220	fOutputList->Add(fV0SemiCentHwAll);
	221
	222	fV0CentTr = new TH1F("fV0CentTr","Centrality percentile with hardware V0 charge thresholds",400,0,100);
	223	fOutputList->Add(fV0CentTr);
	224	fV0CentTrAll = new TH1F("fV0CentTrAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
	225	fOutputList->Add(fV0CentTrAll);
	226	fV0SemiCentTr = new TH1F("fV0SemiCentTr","Centrality percentile with hardware V0 charge thresholds",400,0,100);
	227	fOutputList->Add(fV0SemiCentTr);
	228	fV0SemiCentTrAll = new TH1F("fV0SemiCentTrAll","Centrality percentile with hardware V0 charge thresholds (no event selection)",100,0,100);
	229	fOutputList->Add(fV0SemiCentTrAll);
	230
	231	fV0Percent63 = new TH1F("fV0Percent63","Centrality percentile based on v0 mult",400,0,100);
	232	fOutputList->Add(fV0Percent63);
	233	fV0Percent63All = new TH1F("fV0Percent63All","Centrality percentile based on v0 mult (no event selection)",400,0,100);
	234	fOutputList->Add(fV0Percent63All);
	235
	236	fV0MultAll = new TH2F("fV0MultAll","",250,0,25000,250,0,25000);
	237	fOutputList->Add(fV0MultAll);
	238	fV0Mult63 = new TH2F("fV0Mult63","",250,0,25000,250,0,25000);
	239	fOutputList->Add(fV0Mult63);
	240
e5f60c5a	241	fV0CentVsMult = new TH2F("fV0CentVsMult","",250,0,25000,200,0,100);
	242	fOutputList->Add(fV0CentVsMult);
	243	fV0CentVsCharge = new TH2F("fV0CentVsCharge","",250,0,100000,200,0,100);
	244	fOutputList->Add(fV0CentVsCharge);
	245	fV0CentVsTrCharge = new TH2F("fV0CentVsTrCharge","",250,0,50000,200,0,100);
	246	fOutputList->Add(fV0CentVsTrCharge);
	247	fV0TrChargeVsChargeA = new TH2F("fV0TrChargeVsChargeA","",250,0,75000,250,0,37500);
	248	fOutputList->Add(fV0TrChargeVsChargeA);
	249	fV0TrChargeVsChargeC = new TH2F("fV0TrChargeVsChargeC","",250,0,75000,250,0,37500);
	250	fOutputList->Add(fV0TrChargeVsChargeC);
	251
	252	fAdcPmt = new TH2F("fAdcPmt","",64,-0.5,63.5,400,0,4000);
	253	fOutputList->Add(fAdcPmt);
	254	fAdcPmt0 = new TH2F("fAdcPmt0","",64,-0.5,63.5,400,0,4000);
	255	fOutputList->Add(fAdcPmt0);
	256	fAdcPmt1 = new TH2F("fAdcPmt1","",64,-0.5,63.5,400,0,4000);
	257	fOutputList->Add(fAdcPmt1);
	258	fAdcPmt2 = new TH2F("fAdcPmt2","",64,-0.5,63.5,400,0,4000);
	259	fOutputList->Add(fAdcPmt2);
	260	fAdcPmt3 = new TH2F("fAdcPmt3","",64,-0.5,63.5,400,0,4000);
	261	fOutputList->Add(fAdcPmt3);
	262
c014f07d	263	PostData(1, fOutputList);
	264	}
	265	//________________________________________________________________________
	266	void AliAnaVZEROTrigger::Init()
	267	{
	268	// Nothing here
	269	// ....
	270	}
	271
	272	//________________________________________________________________________
	273	void AliAnaVZEROTrigger::UserExec(Option_t *)
	274	{
	275	// Main loop
	276	// Called for each event
	277
	278
	279	fESD = dynamic_cast<AliESDEvent*>(InputEvent());
	280	if (!fESD) {
	281	printf("ERROR: fESD not available\n");
	282	return;
	283	}
	284
	285	AliESDVZERO* esdV0 = fESD->GetVZEROData();
	286	if (!esdV0) {
	287	Printf("ERROR: esd V0 not available");
	288	return;
	289	}
	290
	291	// Trigger
	292	TString trigStr(fESD->GetFiredTriggerClasses());
	293	if (!trigStr.Contains("-B-")) return;
	294	if (!trigStr.Contains(fMBTrigName.Data())) return;
	295
	296	// Count V0 flags
	297	Int_t nV0A = 0;
	298	Int_t nV0C = 0;
	299	for(Int_t i = 0; i < 32; ++i) {
	300	if (esdV0->GetBBFlag(i)) nV0C++;
	301	if (esdV0->GetBBFlag(i+32)) nV0A++;
	302	}
	303
	304	// Phys sel
	305	Bool_t goodEvent = kTRUE;
e5f60c5a	306	Bool_t isSelected;
e5f60c5a	307	if (fUsePhysSel)
8d939043	308	isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
e5f60c5a	309	else
	310	isSelected = ((esdV0->GetV0ADecision()==1) && (esdV0->GetV0CDecision()==1));
	311
c014f07d	312	if (!isSelected) goodEvent = kFALSE;
c014f07d	313
8d939043	314	const AliESDVertex *primaryVtx = fESD->GetPrimaryVertexSPD();
c014f07d	315	if (!primaryVtx) goodEvent = kFALSE;
	316	if (!primaryVtx->GetStatus()) goodEvent = kFALSE;
	317	Double_t tPrimaryVtxPosition[3];
	318	primaryVtx->GetXYZ(tPrimaryVtxPosition);
e5f60c5a	319	if (goodEvent) {
	320	fZvtx->Fill(tPrimaryVtxPosition[2]);
	321	fXYvtx->Fill(tPrimaryVtxPosition[0],tPrimaryVtxPosition[1]);
	322	}
c014f07d	323	if (TMath::Abs(tPrimaryVtxPosition[2]) > 10.0) goodEvent = kFALSE;
c014f07d	324
e5f60c5a	325	if (goodEvent) {
	326	for(Int_t i = 0; i < 64; ++i) {
	327	fAdcPmt->Fill(i,esdV0->GetAdc(i));
	328	if (tPrimaryVtxPosition[2]>5.0) fAdcPmt0->Fill(i,esdV0->GetAdc(i));
	329	if (tPrimaryVtxPosition[2]>0.0 && tPrimaryVtxPosition[2]<5.0) fAdcPmt1->Fill(i,esdV0->GetAdc(i));
	330	if (tPrimaryVtxPosition[2]>-5.0 && tPrimaryVtxPosition[2]<0.0) fAdcPmt2->Fill(i,esdV0->GetAdc(i));
	331	if (tPrimaryVtxPosition[2]<-5.0) fAdcPmt3->Fill(i,esdV0->GetAdc(i));
	332	}
	333	}
	334
	335	UShort_t chargeA = esdV0->GetTriggerChargeA();
	336	UShort_t chargeC = esdV0->GetTriggerChargeC();
	337
	338	Float_t offChargeA = 0;
	339	Float_t offChargeC = 0;
	340	for(Int_t i = 0; i < 64; ++i) {
	341	if (i < 32) offChargeC += esdV0->GetAdc(i);
	342	else offChargeA += esdV0->GetAdc(i);
	343	}
	344
	345	if (goodEvent) {
	346	fV0TrChargeVsChargeA->Fill(offChargeA,chargeA);
	347	fV0TrChargeVsChargeC->Fill(offChargeC,chargeC);
	348	}
	349
	350	Float_t percentile = 0;
c014f07d	351	AliCentrality *centrality = fESD->GetCentrality();
e5f60c5a	352	// if (centrality->GetQuality()) goodEvent = kFALSE;
	353	percentile = centrality->GetCentralityPercentile("V0M");
	354	// if (percentile < 0) percentile = 0;
	355
	356	if (goodEvent) {
	357	fV0CentVsMult->Fill(esdV0->GetMTotV0A()+esdV0->GetMTotV0C(),percentile);
	358	fV0CentVsCharge->Fill(offChargeA+offChargeC,percentile);
	359	fV0CentVsTrCharge->Fill(chargeA+chargeC,percentile);
	360	}
c014f07d	361
	362	fV0MultAll->Fill(esdV0->GetMTotV0A(),esdV0->GetMTotV0C());
	363	if ((nV0A+nV0C)>=63)
	364	fV0Mult63->Fill(esdV0->GetMTotV0A(),esdV0->GetMTotV0C());
	365
	366	fV0PercentAll->Fill(percentile);
	367	if (goodEvent) fV0Percent->Fill(percentile);
	368
	369	if ((nV0A+nV0C)>=63) {
	370	fV0Percent63All->Fill(percentile);
	371	if (goodEvent) fV0Percent63->Fill(percentile);
	372	}
	373
	374	Float_t multV0 = esdV0->GetMTotV0A()+esdV0->GetMTotV0C();
	375	if (goodEvent) fV0Mult1d->Fill(multV0);
	376
c014f07d	377	fV0Charge2dAll->Fill(chargeA,chargeC);
	378	if (goodEvent) {
	379	fV0Charge2d->Fill(chargeA,chargeC);
	380	fV0Charge2dPercent->Fill(chargeA,chargeC,percentile);
	381	}
	382
	383	for(Int_t j = 0; j < fNThr; ++j) {
	384	Float_t thrA = GetThrA(j);
	385	Float_t thrC = GetThrC(j);
	386	if (chargeA >= thrA && chargeC >= thrC) {
	387	fV0PercentBinsAll[j]->Fill(percentile);
	388	if (goodEvent) fV0PercentBins[j]->Fill(percentile);
	389	}
	390	}
	391
	392	if (chargeA >= fCentCuts[0] && chargeC >= fCentCuts[1]) {
	393	fV0CentAll->Fill(percentile);
	394	if (goodEvent) fV0Cent->Fill(percentile);
	395	}
	396	if (chargeA >= fSemiCentCuts[0] && chargeC >= fSemiCentCuts[1]) {
	397	fV0SemiCentAll->Fill(percentile);
	398	if (goodEvent) fV0SemiCent->Fill(percentile);
	399	}
	400
	401	if (esdV0->GetTriggerBits() & (1<<8)) {
	402	fV0CentHwAll->Fill(percentile);
	403	if (goodEvent) fV0CentHw->Fill(percentile);
	404	}
	405	if (esdV0->GetTriggerBits() & (1<<6)) {
	406	fV0SemiCentHwAll->Fill(percentile);
	407	if (goodEvent) fV0SemiCentHw->Fill(percentile);
	408	}
	409
	410	if (trigStr.Contains("CVHN")) {
b13c9bc7	411	fV0CentTrAll->Fill(percentile);
b13c9bc7	412	if (goodEvent) fV0CentTr->Fill(percentile);
c014f07d	413	}
c014f07d	414	if (trigStr.Contains("CVLN")) {
b13c9bc7	415	fV0SemiCentTrAll->Fill(percentile);
b13c9bc7	416	if (goodEvent) fV0SemiCentTr->Fill(percentile);
c014f07d	417	}
	418
	419	PostData(1, fOutputList);
	420	}
	421
	422	//________________________________________________________________________
	423	void AliAnaVZEROTrigger::Terminate(Option_t *)
	424	{
	425	// Check the output list and store output config file
	426	// Called once at the end of the query
	427
	428	fOutputList = dynamic_cast<TList*> (GetOutputData(1));
	429	if (!fOutputList) {
	430	printf("ERROR: Output list not available\n");
	431	return;
	432	}
	433
	434	printf("\n\n");
	435	FILE *fout=fopen("./trigger.txt","w");
	436	if (!fout) {
	437	printf("Failed to open local result file\n");
	438	return;
	439	}
	440	printf("%.0f %.0f %f %d %.0f %.0f %.0f %.0f\n\n",
	441	GetMinThr(),
	442	GetMaxThr(),
	443	GetRatio(),
	444	GetNThr(),
	445	GetSemiCentCuts(0),GetSemiCentCuts(1),
	446	GetCentCuts(0),GetCentCuts(1));
	447
	448	fprintf(fout,"%.0f %.0f %f %d %.0f %.0f %.0f %.0f\n\n",
	449	GetMinThr(),
	450	GetMaxThr(),
	451	GetRatio(),
	452	GetNThr(),
	453	GetSemiCentCuts(0),GetSemiCentCuts(1),
	454	GetCentCuts(0),GetCentCuts(1));
	455	fclose(fout);
	456	}
	457
	458	//________________________________________________________________________
	459	Float_t AliAnaVZEROTrigger::GetThrA(Int_t j) const
	460	{
	461	// Get the threshold on A side with
	462	// index i
	463	if (j<0 \|\| j>= fNThr) return 0;
	464
	465	return (fMinThr + ((Float_t)j)*(fMaxThr-fMinThr)/((Float_t)fNThr-1.));
	466	}
	467
	468	//________________________________________________________________________
	469	Float_t AliAnaVZEROTrigger::GetThrC(Int_t j) const
	470	{
	471	// Get the threshold on C side with
	472	// index i
	473	if (j<0 \|\| j>= fNThr) return 0;
	474
	475	Float_t thrA = (fMinThr + ((Float_t)j)*(fMaxThr-fMinThr)/((Float_t)fNThr-1.));
	476	return (thrA*fRatio);
	477	}
	478