// Default constructor
fPhysicsSelection = new AliPhysicsSelection;
fPhysicsSelection->SetAnalyzeMC(kTRUE); //For the background correction. This is reset in Init if relevant
+ // fPhysicsSelection->SetUseBXNumbers(kFALSE);
+
AliBackgroundSelection* backgroundSelection = new AliBackgroundSelection("bg","bg");
-
#include <TROOT.h>
#include <TSystem.h>
#include <TInterpreter.h>
#include <TFile.h>
#include <TList.h>
#include <iostream>
+#include <string>
+#include <TCanvas.h>
#include "TH1F.h"
#include "TH2F.h"
+#include "TProfile.h"
#include "AliFMDAnalysisTaskBFCorrelation.h"
#include "AliAnalysisManager.h"
#include "AliESDFMD.h"
//#include "TDatabasePDG.h"
//#include "TParticlePDG.h"
#include "AliESDInputHandler.h"
-ClassImp(AliFMDAnalysisTaskBFCorrelation)
+ClassImp(AliFMDAnalysisTaskBFCorrelation)
AliFMDAnalysisTaskBFCorrelation::AliFMDAnalysisTaskBFCorrelation()
: fDebug(0),
fOutputList(0),
fInputList(0),
+ fInternalList(0),
fVertexString(0x0),
- fStandalone(kTRUE)
+ fStandalone(kTRUE),
+ fEvent(0),
+ fnBinsX(200),
+ fXmin(-6),
+ fXmax(6),
+ fnBinsY(20),
+ fYmin(0),
+ fYmax(2 * TMath::Pi()),
+ c(0),
+ debug0(0),
+ debug1(0)
{
// Default constructor
DefineInput (0, TList::Class());
DefineOutput(0, TList::Class());
}
+
//_____________________________________________________________________
AliFMDAnalysisTaskBFCorrelation::AliFMDAnalysisTaskBFCorrelation(const char* name, Bool_t SE):
AliAnalysisTask(name,name),
- fDebug(0),
- fOutputList(0),
- fInputList(0),
- fVertexString(0x0),
- fStandalone(kTRUE)
+ fDebug(0),
+ fOutputList(0),
+ fInputList(0),
+ fInternalList(0),
+ fVertexString(0x0),
+ fStandalone(kTRUE),
+ fEvent(0),
+ fnBinsX(200),
+ fXmin(-6),
+ fXmax(6),
+ fnBinsY(20),
+ fYmin(0),
+ fYmax(2 * TMath::Pi()),
+ c(0),
+ debug0(0),
+ debug1(0)
{
fStandalone = SE;
if(fStandalone) {
DefineInput (0, TList::Class());
DefineInput(1, TObjString::Class());
DefineOutput(0, TList::Class());
-
}
}
+
//_____________________________________________________________________
void AliFMDAnalysisTaskBFCorrelation::CreateOutputObjects()
{
- //AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
-
+ // Setup the list for storing results, if it does not exist
+
if(!fOutputList) {
fOutputList = new TList();
fOutputList->SetName("BackgroundCorrected");
}
+
+ // Setup the list for temporary storage during the analysis
+
+ if(!fInternalList) {
+ fInternalList = new TList();
+ fInternalList->SetName("InternalBFList");
+ }
+
+ // Set up histograms for analysis and storage. 4 different binnings
- TH1F* test = new TH1F("test","test",10,0,10);
- fOutputList->Add(test);
+ for (Int_t i = 1; i <= 4; i++) {
+
+ // Temporary histograms for storing hist per eta summed over phi
+
+ TH1D *hESDMult = new TH1D(Form("hESDMult_binning%d", i),
+ Form("Multiplicity per event vs eta ESD (%d bins)", 30*i),
+ 18*i, -9, 9);
+ TH1D *hMCMult = new TH1D(Form("hMCMult_binning%d", i),
+ Form("Multiplicity per event vs eta MC-truth (%d bins)", 30*i),
+ 18*i, -9, 9);
+ fInternalList->Add(hESDMult);
+ fInternalList->Add(hMCMult);
+
+
+ // Histograms for storing the 5 values to calculate b (ESD & MC)
+
+ TH1D *pnFnB = new TH1D(Form("pnFnB_binning%d", i),
+ Form("Correlation ESD (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnF2 = new TH1D(Form("pnF2_binning%d", i),
+ Form("Fwd^2 ESD (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnB2 = new TH1D(Form("pnB2_binning%d", i),
+ Form("Bwd^2 ESD (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnF = new TH1D(Form("pnF_binning%d", i),
+ Form("Fwd ESD (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnB = new TH1D(Form("pnB_binning%d", i),
+ Form("Bwd ESD (%d bins)", i*9),
+ 9*i, 0, 9);
+
+ TH1D *pnFnB_MC = new TH1D(Form("pnFnB_MC_binning%d", i),
+ Form("Correlation MC (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnF2_MC = new TH1D(Form("pnF2_MC_binning%d", i),
+ Form("Fwd^2 MC (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnB2_MC = new TH1D(Form("pnB2_MC_binning%d", i),
+ Form("Bwd^2 MC (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnF_MC = new TH1D(Form("pnF_MC_binning%d", i),
+ Form("Fwd MC (%d bins)", i*9),
+ 9*i, 0, 9);
+ TH1D *pnB_MC = new TH1D(Form("pnB_MC_binning%d", i),
+ Form("Backward MC (%d bins)" ,i*9),
+ 9*i, 0, 9);
+ fOutputList->Add(pnFnB);
+ fOutputList->Add(pnF2);
+ fOutputList->Add(pnB2);
+ fOutputList->Add(pnF);
+ fOutputList->Add(pnB);
+
+ fOutputList->Add(pnFnB_MC);
+ fOutputList->Add(pnF2_MC);
+ fOutputList->Add(pnB2_MC);
+ fOutputList->Add(pnF_MC);
+ fOutputList->Add(pnB_MC);
+
+ // Debugging histograms : "dNdEta"
+
+ // TH1D *hESDMultvsEta = new TH1D(Form("hESDMultvsEta_binning%d" ,i),
+ // Form("Mult vs Eta (%d bins, ESD)", 18*i),
+ // 18*i, -9, 9);
+ TH1D *hESDMultvsEta = new TH1D(Form("hESDMultvsEta_binning%d" ,i),
+ Form("Mult vs Eta (%d bins, ESD)", i*18),
+ 200, -6, 6);
+ TH1D *hMCMultvsEta = new TH1D(Form("hMCMultvsEta_binning%d", i),
+ Form("Mult vs Eta (%d bins, MC)", 18*i),
+ 200, -6, 6);
+ //TH1D *hMCMultvsEta = new TH1D(Form("hMCMultvsEta_binning%d", i),
+ // Form("Mult vs Eta (%d bins, MC)", 18*i),
+ // 18*i, -9, 9);
+ fOutputList->Add(hESDMultvsEta);
+ fOutputList->Add(hMCMultvsEta);
+
+ TH1D *hESDMultvsEtaC = new TH1D(Form("hESDMultvsEtaC_binning%d" ,i),
+ Form("Mult vs Eta C (%d bins, ESD)", 18*i),
+ 18*i, -9, 9);
+ TH1D *hMCMultvsEtaC = new TH1D(Form("hMCMultvsEtaC_binning%d", i),
+ Form("Mult vs Eta C (%d bins, MC)", 18*i),
+ 18*i, -9, 9);
+ fOutputList->Add(hESDMultvsEtaC);
+ fOutputList->Add(hMCMultvsEtaC);
+
+ // Debugging histograms : Different distributions
+
+ TH2D *hESDBwdDist = new TH2D(Form("hESDBwdDist_binning%d", i),
+ Form("Distribution of Bwd values (%d bins)", i*9),
+ 9*i, 0, 9, 20, 0, 20);
+
+ TH2D *hESDBwd2Dist = new TH2D(Form("hESDBwd2Dist_binning%d", i),
+ Form("Distribution of Bwd^2 values (%d bins)", i*9),
+ 9*i, 0, 9, 400, 0, 400);
+
+ TH2D *hMCBwdDist = new TH2D(Form("hMCBwdDist_binning%d", i),
+ Form("Distribution of Bwd values (%d bins)", i*9),
+ 9*i, 0, 9, 20, 0, 20);
+
+ TH2D *hMCBwd2Dist = new TH2D(Form("hMCBwd2Dist_binning%d", i),
+ Form("Distribution of Bwd^2 values (%d bins)", i*9),
+ 9*i, 0, 9, 400, 0, 400);
+ fOutputList->Add(hESDBwdDist);
+ fOutputList->Add(hESDBwd2Dist);
+ fOutputList->Add(hMCBwdDist);
+ fOutputList->Add(hMCBwd2Dist);
+ }
}
+
//_____________________________________________________________________
void AliFMDAnalysisTaskBFCorrelation::ConnectInputData(Option_t */*option*/)
{
fVertexString = (TObjString*)GetInputData(1);
}
}
+
//_____________________________________________________________________
void AliFMDAnalysisTaskBFCorrelation::Exec(Option_t */*option*/)
{
+ fEvent++;
+ if (fEvent % 1000 == 0)
+ std::cout << "Event # " << fEvent << std::endl;
+
AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
fVertexString = (TObjString*)fInputList->At(0);
-
Int_t vtxbin = fVertexString->GetString().Atoi();
- //for(UShort_t det=1;det<=3;det++) {
- // Int_t nRings = (det==1 ? 1 : 2);
- // for (UShort_t ir = 0; ir < nRings; ir++) {
- // Char_t ringChar = (ir == 0 ? 'I' : 'O');
- //TH2F* hMultTotal = (TH2F*)fOutputList->FindObject(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
- // TH2F* hMultTotalTrVtx = (TH2F*)fOutputList->FindObject(Form("dNdetaTrVtx_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
-
+ if(vtxbin !=4) return;
+ SetBounds();
+ CountESDHits();
+ CalculateValues("ESD");
- TH1F* test = (TH1F*)fOutputList->FindObject("test");
- test->Fill(vtxbin);
if(pars->GetProcessPrimary())
ProcessPrimary();
if(fStandalone) {
PostData(0, fOutputList);
}
+}
+//_____________________________________________________________________
+void AliFMDAnalysisTaskBFCorrelation::CountESDHits() {
+
+ TH2F *dNdEtadPhiHist = (TH2F*)fInputList->FindObject("dNdetadphiHistogramSPDTrVtx");
+
+ for (Int_t i = 1; i<=4; i++) {
+
+ TH1D *hESDMult = (TH1D*)fInternalList->FindObject(Form("hESDMult_binning%d", i));
+ hESDMult->Reset();
+ }
+
+ Int_t etamax = dNdEtadPhiHist->GetNbinsX();
+ Int_t phimax = dNdEtadPhiHist->GetNbinsY();
+
+ TH1D *hprod = dNdEtadPhiHist->ProjectionX();
+
+ TH1D *hESDMultvsEta = (TH1D*)fOutputList->FindObject(Form("hESDMultvsEta_binning%d" ,4));
+ hESDMultvsEta->Add(hprod);
+
+ for (Int_t etabin = 1; etabin <= etamax; etabin++) {
+ Float_t val = 0;
+ Float_t eta = dNdEtadPhiHist->GetXaxis()->GetBinCenter(etabin);
+
+
+ for (Int_t i = 1; i <= 4; i++) {
+ TH1D *hESDMult = (TH1D*)fInternalList->FindObject(Form("hESDMult_binning%d", i));
+ hESDMult->Fill(eta, hprod->GetBinContent(etabin));
+ }
+ }
+
+ //hESDMult->Add(hprod);
}
//_____________________________________________________________________
-void AliFMDAnalysisTaskBFCorrelation::Terminate(Option_t */*option*/) {
+void AliFMDAnalysisTaskBFCorrelation::CalculateValues(TString type) {
+
+ type.ToUpper();
+
+ const char *stype = type.Data();
+
+ TString sinput;
+ TString soutput;
+ TString soutputC;
+ TString snFnB;
+ TString snF2;
+ TString snB2;
+ TString snF;
+ TString snB;
+
+ if (type == "ESD") {
+ snFnB = "pnFnB";
+ snF2 = "pnF2";
+ snB2 = "pnB2";
+ snF = "pnF";
+ snB = "pnB";
+ } else {
+ snFnB = "pnFnB_MC";
+ snF2 = "pnF2_MC";
+ snB2 = "pnB2_MC";
+ snF = "pnF_MC";
+ snB = "pnB_MC";
+ }
+
+ for (Int_t i = 1; i <= 4; i++) {
+
+ sinput = Form("h%sMult_binning%d", stype, i);
+ soutput = Form("h%sMultvsEta_binning%d", stype, i);
+ soutputC = Form("h%sMultvsEtaC_binning%d", stype, i);
+
+ TH1D *input = (TH1D*)fInternalList->FindObject(sinput);
+ TH1D *output = (TH1D*)fOutputList->FindObject(soutput);
+ TH1D *outputC = (TH1D*)fOutputList->FindObject(soutputC);
+
+ //output->Add(input);
+
+ for (Int_t bin = 1; bin <= input->GetNbinsX()/2; bin++) {
+
+ Double_t eta = input->GetBinCenter(bin);
+ Double_t bwd = input->GetBinContent(bin);
+ Double_t bwd2 = TMath::Power(bwd, 2);
+ Double_t fwd = input->GetBinContent(input->FindBin(TMath::Abs(eta)));
+ Double_t fwd2 = TMath::Power(fwd, 2);
+ Double_t fwdbwd = fwd*bwd;
+
+ snFnB += Form("_binning%d", i);
+ snF2 += Form("_binning%d", i);
+ snB2 += Form("_binning%d", i);
+ snF += Form("_binning%d", i);
+ snB += Form("_binning%d", i);
+
+ TH1D *nFnB = static_cast<TH1D*>(fOutputList->FindObject(snFnB));
+ TH1D *nF2 = static_cast<TH1D*>(fOutputList->FindObject(snF2));
+ TH1D *nB2 = static_cast<TH1D*>(fOutputList->FindObject(snB2));
+ TH1D *nF = static_cast<TH1D*>(fOutputList->FindObject(snF));
+ TH1D *nB = static_cast<TH1D*>(fOutputList->FindObject(snB));
+
+ nFnB->Fill(TMath::Abs(eta), fwdbwd);
+ nF2->Fill(TMath::Abs(eta), fwd2);
+ nB2->Fill(TMath::Abs(eta), bwd2);
+ nF->Fill(TMath::Abs(eta), fwd);
+ nB->Fill(TMath::Abs(eta), bwd);
+
+ outputC->Fill(TMath::Abs(eta), fwd);
+ outputC->Fill(eta, bwd);
+
+ if (type == "ESD") {
+
+ snFnB.Remove(5,9);
+ snF2.Remove(4,9);
+ snB2.Remove(4,9);
+ snF.Remove(3,9);
+ snB.Remove(3,9);
+ TH2D *hESDBwdDist = (TH2D*)fOutputList->FindObject(Form("hESDBwdDist_binning%d", i));
+ TH2D *hESDBwd2Dist = (TH2D*)fOutputList->FindObject(Form("hESDBwd2Dist_binning%d", i));
+ hESDBwdDist->Fill(TMath::Abs(eta), bwd);
+ hESDBwd2Dist->Fill(TMath::Abs(eta), bwd2);
+
+ } else {
+
+ snFnB.Remove(8,9);
+ snF2.Remove(7,9);
+ snB2.Remove(7,9);
+ snF.Remove(6,9);
+ snB.Remove(6,9);
+ TH2D *hMCBwdDist = (TH2D*)fOutputList->FindObject(Form("hMCBwdDist_binning%d", i));
+ TH2D *hMCBwd2Dist = (TH2D*)fOutputList->FindObject(Form("hMCBwd2Dist_binning%d", i));
+ hMCBwdDist->Fill(TMath::Abs(eta), bwd);
+ hMCBwd2Dist->Fill(TMath::Abs(eta), bwd2);
+ }
+ }
+ }
+}
+
+//_____________________________________________________________________
+void AliFMDAnalysisTaskBFCorrelation::SetBounds() {
+
+ TH2D *hTemp = (TH2D*)fInputList->FindObject("multTrVtx_FMD1I_vtxbin0");
+ fnBinsX = hTemp->GetNbinsX();
+ fnBinsY = hTemp->GetNbinsY();
+ fXmin = hTemp->GetXaxis()->GetXmin();
+ fXmax = hTemp->GetXaxis()->GetXmax();
+ fYmin = hTemp->GetYaxis()->GetXmin();
+ fYmax = hTemp->GetYaxis()->GetXmax();
+}
+
+//_____________________________________________________________________
+void AliFMDAnalysisTaskBFCorrelation::Terminate(Option_t */*option*/) {
}
//_____________________________________________________________________
void AliFMDAnalysisTaskBFCorrelation::ProcessPrimary() {
- /*
+
AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
AliMCEvent* mcEvent = eventHandler->MCEvent();
if(!mcEvent)
return;
-
AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
AliHeader* header = mcEvent->Header();
AliGenEventHeader* genHeader = header->GenEventHeader();
+ AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
+
+ if (!pythiaGenHeader) {
+ std::cout<<" no pythia header!"<<std::endl;
+ return;
+ }
+
+ Int_t pythiaType = pythiaGenHeader->ProcessType();
+
+ if(pythiaType==92||pythiaType==93){
+ std::cout<<"single diffractive"<<std::endl;
+ return;
+ }
+ if(pythiaType==94){
+ std::cout<<"double diffractive"<<std::endl;
+ return;
+ }
+
+
TArrayF vertex;
genHeader->PrimaryVertex(vertex);
if(TMath::Abs(vertex.At(2)) > pars->GetVtxCutZ())
return;
- Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
- Double_t vertexBinDouble = (vertex.At(2) + pars->GetVtxCutZ()) / delta;
- Int_t vertexBin = (Int_t)vertexBinDouble;
+ //Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
+ //Double_t vertexBinDouble = (vertex.At(2) + pars->GetVtxCutZ()) / delta;
+ //Int_t vertexBin = (Int_t)vertexBinDouble;
- Bool_t firstTrack = kTRUE;
+ //Bool_t firstTrack = kTRUE;
// we loop over the primaries only unless we need the hits (diagnostics running slowly)
Int_t nTracks = stack->GetNprimary();
- if(pars->GetProcessHits())
- nTracks = stack->GetNtrack();
+ // if(pars->GetProcessHits())
+ // nTracks = stack->GetNtrack();
+ TH1D *hMCMultvsEta = (TH1D*)fOutputList->FindObject(Form("hMCMultvsEta_binning%d" ,4));
- for(Int_t i = 0 ;i<nTracks;i++) {
- particle = (AliMCParticle*) mcEvent->GetTrack(i);
- if(!particle)
+ for (Int_t i = 1; i <= 4; i++) {
+
+ TH1D *hMCMult = (TH1D*)fInternalList->FindObject(Form("hMCMult_binning%d", i));
+ hMCMult->Reset();
+
+ for(Int_t ii = 0 ;ii<nTracks;ii++) {
+ particle = (AliMCParticle*) mcEvent->GetTrack(ii);
+ if(!particle)
continue;
-
- if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
- hPrimary->Fill(particle->Eta());
-
- TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin));
- hPrimVtxBin->Fill(particle->Eta());
- if(firstTrack) {
- fNMCevents.Fill(vertexBin);
- firstTrack = kFALSE;
+ if(stack->IsPhysicalPrimary(ii) && particle->Charge() != 0) {
+ if (i == 1)
+ hPrimary->Fill(particle->Eta());
+ if(i==4)
+ hMCMultvsEta->Fill(particle->Eta());
+ hMCMult->Fill(particle->Eta());
}
}
}
- */
+ CalculateValues("MC");
}
//_____________________________________________________________________
//
#include "TObjString.h"
#include "TArrayI.h"
#include "TH1I.h"
+#include "TH2.h"
#include "AliMCEvent.h"
#include "AliFMDFloatMap.h"
+#include "TCanvas.h"
/**
* @ingroup FMD_ana
fDebug(o.fDebug),
fOutputList(0),
fInputList(0),
+ fInternalList(0),
fVertexString(o.fVertexString),
- fStandalone(o.fStandalone)
- {}
+ fStandalone(o.fStandalone),
+ fEvent(0),
+ fnBinsX(200),
+ fXmin(-6),
+ fXmax(6),
+ fnBinsY(20),
+ fYmin(2),
+ fYmax(2 * TMath::Pi()),
+ c(0),
+ debug0(0),
+ debug1(0)
+ {}
AliFMDAnalysisTaskBFCorrelation& operator=(const AliFMDAnalysisTaskBFCorrelation&) { return *this; }
// Implementation of interface methods
void SetInputList(TList* inputList) {fInputList = inputList;}
void SetInputVertex(TObjString* vtxString) {fVertexString = vtxString;}
void SetOutputList(TList* outputList) {fOutputList = outputList;}
+ void CountESDHits();
+ void CalculateValues(TString type);
+ void SetBounds();
+ void MergeEvent(TH2D *hMultTrVtxFull);
void ProcessPrimary();
TList* GetOutputList() {return fOutputList;}
private:
+
Int_t fDebug; // Debug flag
TList* fOutputList;
TList* fInputList;
+ TList* fInternalList;
TObjString* fVertexString;
Bool_t fStandalone;
+
+ Int_t fEvent;
+ Int_t fnBinsX;
+ Float_t fXmin;
+ Float_t fXmax;
+ Int_t fnBinsY;
+ Float_t fYmin;
+ Float_t fYmax;
+
+ TCanvas *c;
+
+ Double_t debug0;
+ Double_t debug1;
+
ClassDef(AliFMDAnalysisTaskBFCorrelation, 0); // Analysis task for FMD analysis
};
+++ /dev/null
-#include "AliFMDAnalysisTaskGenerateBackground.h"
-#include "AliESDEvent.h"
-#include "iostream"
-#include "AliESDFMD.h"
-#include "TH2F.h"
-#include "AliTrackReference.h"
-#include "AliStack.h"
-#include "AliFMDAnaParameters.h"
-#include "AliFMDStripIndex.h"
-#include "AliStack.h"
-#include "AliMCParticle.h"
-#include "AliMCEvent.h"
-//#include "AliFMDGeometry.h"
-#include "TArray.h"
-#include "AliGenEventHeader.h"
-#include "AliHeader.h"
-#include "AliFMDAnaCalibBackgroundCorrection.h"
-//#include "AliCDBManager.h"
-//#include "AliCDBId.h"
-//#include "AliCDBMetaData.h"
-#include "TSystem.h"
-#include "TROOT.h"
-#include "TAxis.h"
-ClassImp(AliFMDAnalysisTaskGenerateBackground)
-
-//_____________________________________________________________________
-AliFMDAnalysisTaskGenerateBackground::AliFMDAnalysisTaskGenerateBackground():
-AliAnalysisTaskSE(),
- fListOfHits(),
- fListOfPrimaries(),
- fListOfCorrection(),
- fVertexBins(),
- fLastTrackByStrip(0),
- fHitsByStrip(0),
- fZvtxCut(10),
- fNvtxBins(10),
- fNbinsEta(200),
- fBackground(0)
-{
- // Default constructor
-}
-//_____________________________________________________________________
-AliFMDAnalysisTaskGenerateBackground::AliFMDAnalysisTaskGenerateBackground(const char* name):
- AliAnalysisTaskSE(name),
- fListOfHits(),
- fListOfPrimaries(),
- fListOfCorrection(),
- fVertexBins(),
- fLastTrackByStrip(0),
- fHitsByStrip(0),
- fZvtxCut(10),
- fNvtxBins(10),
- fNbinsEta(200),
- fBackground(0)
-{
-
- DefineOutput(1, TList::Class());
- DefineOutput(2, TList::Class());
- DefineOutput(3, TH1F::Class());
- DefineOutput(4, TList::Class());
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::UserCreateOutputObjects()
-{
-// Create the output containers
-//
-
- std::cout<<"Creating output objects"<<std::endl;
- for(Int_t iring = 0; iring<2;iring++) {
- Char_t ringChar = (iring == 0 ? 'I' : 'O');
- Int_t nSec = (iring == 1 ? 40 : 20);
- for(Int_t v=0; v<fNvtxBins;v++) {
-
- TH2F* hPrimary = new TH2F(Form("hPrimary_FMD_%c_vtx%d",ringChar,v),
- Form("hPrimary_FMD_%c_vtx%d",ringChar,v),
- fNbinsEta, -6,6, nSec, 0,2*TMath::Pi());
- hPrimary->Sumw2();
- fListOfPrimaries.Add(hPrimary);
- }
- }
-
-
- for(Int_t det =1; det<=3;det++) {
- Int_t nRings = (det==1 ? 1 : 2);
- for(Int_t ring = 0;ring<nRings;ring++) {
- Int_t nSec = (ring == 1 ? 40 : 20);
- Char_t ringChar = (ring == 0 ? 'I' : 'O');
- TH1F* doubleHits = new TH1F(Form("DoubleHits_FMD%d%c",det,ringChar),
- Form("DoubleHits_FMD%d%c",det,ringChar),fNbinsEta, -6,6);
- TH1F* allHits = new TH1F(Form("allHits_FMD%d%c",det,ringChar),
- Form("allHits_FMD%d%c",det,ringChar), fNbinsEta, -6,6);
-
- doubleHits->Sumw2();
- allHits->Sumw2();
- fListOfHits.Add(allHits);
- fListOfHits.Add(doubleHits);
-
- for(Int_t v=0; v<fNvtxBins;v++) {
- TH2F* hHits = new TH2F(Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
- Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
- fNbinsEta, -6,6, nSec, 0, 2*TMath::Pi());
- hHits->Sumw2();
- fListOfHits.Add(hHits);
-
- }
- }
- }
-
- fVertexBins.SetName("VertexBins");
- fVertexBins.GetXaxis()->Set(fNvtxBins,-1*fZvtxCut,fZvtxCut);
-
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::Init()
-{
- fLastTrackByStrip.Reset(-1);
-
-
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::UserExec(Option_t */*option*/)
-{
-
- fLastTrackByStrip.Reset(-1);
- fHitsByStrip.Reset(0);
- AliMCEvent* mcevent = MCEvent();
- AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
-
- AliMCParticle* particle = 0;
- AliStack* stack = mcevent->Stack();
-
- UShort_t det,sec,strip;
- Char_t ring;
-
- Int_t nTracks = mcevent->GetNumberOfTracks();
- AliHeader* header = mcevent->Header();
- AliGenEventHeader* genHeader = header->GenEventHeader();
-
- TArrayF vertex;
- genHeader->PrimaryVertex(vertex);
-
- if(TMath::Abs(vertex.At(2)) > fZvtxCut)
- return;
- for(Int_t i = 0 ;i<nTracks;i++) {
- particle = mcevent->GetTrack(i);
-
- if(!particle)
- continue;
-
- Double_t delta = 2*fZvtxCut/fNvtxBins;
- Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta;
- Int_t vertexBin = (Int_t)vertexBinDouble;
-
- if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
-
-
- TH2F* hPrimaryInner = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",'I',vertexBin));
- TH2F* hPrimaryOuter = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",'O',vertexBin));
- hPrimaryInner->Fill(particle->Eta(),particle->Phi());
- hPrimaryOuter->Fill(particle->Eta(),particle->Phi());
- }
-
- for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) {
-
- AliTrackReference* ref = particle->GetTrackReference(j);
-
- if(ref->DetectorId() != AliTrackReference::kFMD)
- continue;
- AliFMDStripIndex::Unpack(ref->UserId(),det,ring,sec,strip);
- Float_t thisStripTrack = fLastTrackByStrip(det,ring,sec,strip);
- if(particle->Charge() != 0 && i != thisStripTrack ) {
- //Double_t x,y,z;
- //AliFMDGeometry* fmdgeo = AliFMDGeometry::Instance();
- //fmdgeo->Detector2XYZ(det,ring,sec,strip,x,y,z);
- Float_t phi = pars->GetPhiFromSector(det,ring,sec);
- Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));
- //Float_t phi = TMath::ATan2(y,x);
- //if(phi<0) phi = phi+2*TMath::Pi();
- // Float_t r = TMath::Sqrt(TMath::Power(x,2)+TMath::Power(y,2));
- //Float_t theta = TMath::ATan2(r,z-vertex.At(2));
- //Float_t eta = -1*TMath::Log(TMath::Tan(0.5*theta));
- TH2F* hHits = (TH2F*)fListOfHits.FindObject(Form("hHits_FMD%d%c_vtx%d", det,ring,vertexBin));
- hHits->Fill(eta,phi);
- Float_t nstrips = (ring =='O' ? 256 : 512);
- fHitsByStrip(det,ring,sec,strip) +=1;
- TH1F* allHits = (TH1F*)fListOfHits.FindObject(Form("allHits_FMD%d%c",det,ring));
- TH1F* doubleHits = (TH1F*)fListOfHits.FindObject(Form("DoubleHits_FMD%d%c",det,ring));
-
- if(fHitsByStrip(det,ring,sec,strip) == 1)
- allHits->Fill(eta);
-
- doubleHits->Fill(eta);
- /*if(fHitsByStrip.operator()(det,ring,sec,strip) == 2){
- TH1F* doubleHits = (TH1F*)fListOfHits.FindObject("DoubleHits");
- doubleHits->Fill(eta,2);
- }*/
- //if(fHitsByStrip.operator()(det,ring,sec,strip) > 1){
- // doubleHits->Fill(eta);
- // }
-
- fLastTrackByStrip(det,ring,sec,strip) = (Float_t)i;
- if(strip >0)
- fLastTrackByStrip(det,ring,sec,strip-1) = (Float_t)i;
- if(strip < (nstrips - 1))
- fLastTrackByStrip(det,ring,sec,strip+1) = (Float_t)i;
- }
- }
-
- }
-
-
- PostData(1, &fListOfHits);
- PostData(2, &fListOfPrimaries);
- PostData(3, &fVertexBins);
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::Terminate(Option_t */*option*/)
-{
- /* TH1F* allHits = (TH1F*)fListOfHits.FindObject("allHits");
- TH1F* doubleHits = (TH1F*)fListOfHits.FindObject("DoubleHits");
-
- doubleHits->Divide(allHits);
- GenerateCorrection();
- PostData(1, &fListOfHits);
- PostData(4, &fListOfCorrection);*/
-
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::GenerateCorrection() {
-
- fBackground = new AliFMDAnaCalibBackgroundCorrection();
-
- for(Int_t det= 1; det <=3; det++) {
- Int_t nRings = (det==1 ? 1 : 2);
-
- for(Int_t iring = 0; iring<nRings; iring++) {
- Char_t ring = (iring == 0 ? 'I' : 'O');
- TH1F* allHits = (TH1F*)fListOfHits.FindObject(Form("allHits_FMD%d%c",det,ring));
- TH1F* doubleHits = (TH1F*)fListOfHits.FindObject(Form("DoubleHits_FMD%d%c",det,ring));
- fBackground->SetDoubleHitCorrection(det,ring,doubleHits);
- doubleHits->Divide(allHits);
- for(Int_t vertexBin=0;vertexBin<fNvtxBins ;vertexBin++) {
- TH2F* hHits = (TH2F*)fListOfHits.FindObject(Form("hHits_FMD%d%c_vtx%d", det,ring,vertexBin));
- TH2F* hPrimary = (TH2F*)fListOfPrimaries.FindObject( Form("hPrimary_FMD_%c_vtx%d",ring,vertexBin));
- TH2F* hCorrection = (TH2F*)hHits->Clone(Form("FMD%d%c_vtxbin_%d_correction",det,ring,vertexBin));
- hCorrection->Divide(hPrimary);
- hCorrection->SetTitle(hCorrection->GetName());
- fListOfCorrection.Add(hCorrection);
- fBackground->SetBgCorrection(det,ring,vertexBin,hCorrection);
-
-
- }
-
- }
- }
- TAxis refAxis(fNvtxBins,-1*fZvtxCut,fZvtxCut);
- fBackground->SetRefAxis(&refAxis);
-
-}
-//_____________________________________________________________________
-void AliFMDAnalysisTaskGenerateBackground::ReadFromFile(const Char_t* filename, Bool_t storeInOCDB, Int_t /*runNo*/) {
-
- TFile infile(filename);
- TH1F* hVertex = (TH1F*)infile.Get("VertexBins");
- fZvtxCut = hVertex->GetXaxis()->GetXmax();
- fNvtxBins = hVertex->GetXaxis()->GetNbins();
- fVertexBins.SetName("VertexBins");
- fVertexBins.GetXaxis()->Set(fNvtxBins,-1*fZvtxCut,fZvtxCut);
-
- TList* listOfHits = (TList*)infile.Get("Hits");
- TList* listOfPrim = (TList*)infile.Get("Primaries");
-
- for(Int_t det =1; det<=3;det++)
- {
- Int_t nRings = (det==1 ? 1 : 2);
- for(Int_t ring = 0;ring<nRings;ring++)
- {
- Char_t ringChar = (ring == 0 ? 'I' : 'O');
- TH1F* allHits = (TH1F*)listOfHits->FindObject(Form("allHits_FMD%d%c",det,ringChar));
- TH1F* doubleHits = (TH1F*)listOfHits->FindObject(Form("DoubleHits_FMD%d%c",det,ringChar));
- fListOfHits.Add(allHits);
- fListOfHits.Add(doubleHits);
- for(Int_t v=0; v<fNvtxBins;v++)
- {
-
- TH2F* hHits = (TH2F*)listOfHits->FindObject(Form("hHits_FMD%d%c_vtx%d", det,ringChar,v));
- fListOfHits.Add(hHits);
- }
- }
- }
- for(Int_t iring = 0; iring<2;iring++) {
- Char_t ringChar = (iring == 0 ? 'I' : 'O');
- for(Int_t v=0; v<fNvtxBins;v++) {
-
- TH2F* hPrimary = (TH2F*)listOfPrim->FindObject( Form("hPrimary_FMD_%c_vtx%d",ringChar,v));
- fListOfPrimaries.Add(hPrimary);
-
- }
- }
- GenerateCorrection();
-
- TFile fout("backgroundFromFile.root","recreate");
- fListOfHits.Write();
- fListOfPrimaries.Write();
- fListOfCorrection.Write();
- fVertexBins.Write();
-
- fout.Close();
-
- if(storeInOCDB) {
- TFile fcalib("$ALICE_ROOT/FMD/Correction/Background/background.root","RECREATE");
- fBackground->Write(AliFMDAnaParameters::GetBackgroundID());
- fcalib.Close();
- /* AliCDBManager* cdb = AliCDBManager::Instance();
- cdb->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
- AliCDBId id(AliFMDAnaParameters::GetBackgroundPath(),runNo,999999999);
-
- AliCDBMetaData* meta = new AliCDBMetaData;
- meta->SetResponsible(gSystem->GetUserInfo()->fRealName.Data());
- meta->SetAliRootVersion(gROOT->GetVersion());
- meta->SetBeamPeriod(1);
- meta->SetComment("Background Correction for FMD");
- meta->SetProperty("key1", fBackground );
- cdb->Put(fBackground, id, meta);
- */
- }
-
-}
-//_____________________________________________________________________
-//
-// EOF
-//
+++ /dev/null
-#ifndef ALIFMDANALYSISTASKGENERATEBACKGROUND_H
-#define ALIFMDANALYSISTASKGENERATEBACKGROUND_H
-
-#include "AliAnalysisTaskSE.h"
-#include "TList.h"
-#include "AliFMDFloatMap.h"
-#include "TH1F.h"
-
-/**
- * Make a background distribution from simulated data
- * @ingroup FMD_ana
- *
- *
- */
-
-class AliFMDAnaCalibBackgroundCorrection;
-
-class AliFMDAnalysisTaskGenerateBackground : public AliAnalysisTaskSE
-{
- public:
- AliFMDAnalysisTaskGenerateBackground();
- AliFMDAnalysisTaskGenerateBackground(const char* name);
- ~AliFMDAnalysisTaskGenerateBackground() {;}
- AliFMDAnalysisTaskGenerateBackground(const AliFMDAnalysisTaskGenerateBackground& o) : AliAnalysisTaskSE(),
- fListOfHits(),
- fListOfPrimaries(),
- fListOfCorrection(),
- fVertexBins(o.fVertexBins),
- fLastTrackByStrip(o.fLastTrackByStrip),
- fHitsByStrip(o.fHitsByStrip),
- fZvtxCut(o.fZvtxCut),
- fNvtxBins(o.fNvtxBins),
- fNbinsEta(o.fNbinsEta),
- fBackground(o.fBackground)
- {}
- AliFMDAnalysisTaskGenerateBackground& operator=(const AliFMDAnalysisTaskGenerateBackground&) { return *this; }
-
- virtual void Init();
- virtual void UserCreateOutputObjects();
- virtual void UserExec(Option_t* /*option*/);
- void Terminate(Option_t */*option*/);
- void SetZvtxCut(Float_t vtxcut) {fZvtxCut = vtxcut;}
- void SetNvtxBins(Int_t nvtxbins) {fNvtxBins = nvtxbins;}
- void SetNbinsEta(Int_t netabins) {fNbinsEta = netabins;}
- void ReadFromFile(const Char_t* filename = "background.root", Bool_t storeInOCDB = kFALSE, Int_t runNo=0);
- private:
-
- void GenerateCorrection();
-
- TList fListOfHits;
- TList fListOfPrimaries;
- TList fListOfCorrection;
- TH1F fVertexBins;
- AliFMDFloatMap fLastTrackByStrip;
- AliFMDFloatMap fHitsByStrip;
- Float_t fZvtxCut;
- Int_t fNvtxBins;
- Int_t fNbinsEta;
- AliFMDAnaCalibBackgroundCorrection* fBackground;
- ClassDef(AliFMDAnalysisTaskGenerateBackground, 1);
-
-};
-#endif
-// Local Variables:
-// mode: C++
-// End:
if(!vtxStatus)
vtxFound = kFALSE;
- if(TMath::Abs(vertex.At(2)) > fZvtxCut) {
- vtxFound = kFALSE;
+ //if(TMath::Abs(vertex.At(2)) > fZvtxCut) {
+ // vtxFound = kFALSE;
- }
+ //}
Bool_t isTriggered = pars->IsEventTriggered(AliFMDAnaParameters::kMB1);
Bool_t isTriggeredNSD = pars->IsEventTriggered(AliFMDAnaParameters::kNSD);
if(vtxFound && isTriggered) hEventsSelected->Fill(vertexBin);
// if(!vtxFound || !isTriggered) return;
+ if(TMath::Abs(vertex.At(2)) > fZvtxCut) {
+ return;
+ }
+
for(Int_t i = 0 ;i<nTracks;i++) {
particle = (AliMCParticle*) mcevent->GetTrack(i);
nNonZero++;
}
Int_t nBinsOld = fNbinsToCut;
- if(det == 1 && ringChar =='I') {
- fNbinsToCut = 0;
- }
+ //if(det == 1 && ringChar =='I') {
+ // fNbinsToCut = 0;
+ // }
TH1F* hRingMult = (TH1F*)fMultList[what]->FindObject(Form("hRingMult_FMD%d%c_%s",det,ringChar,fAnalysisNames[what]));
for(Int_t i=1; i<=hRingMult->GetNbinsX(); i++) {
else if(pars->GetEnergy() == AliFMDAnaParameters::k900)
fpyt = TFile::Open("/home/canute/ALICE/FMDanalysis/FirstAnalysis/pythia_study/pythiahists900.root","READ");
+
+
if(realdata ) {
if(fpyt) {
hPythiaMB = (TH1F*)fpyt->Get("hPythiaMB");
--- /dev/null
+void SubmitFMDCorrections(const Char_t* filename, Bool_t store, Float_t energy, Int_t trigger, Float_t mag, Int_t collsystem) {
+
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libPWG0base");
+ gSystem->Load("libPWG2forward");
+
+ AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
+
+ if(energy == 900)
+ pars->SetEnergy(AliFMDAnaParameters::k900);
+ else if(energy == 7000)
+ pars->SetEnergy(AliFMDAnaParameters::k7000);
+ else if(energy == 10000)
+ pars->SetEnergy(AliFMDAnaParameters::k10000);
+ else if(energy == 14000)
+ pars->SetEnergy(AliFMDAnaParameters::k14000);
+ else if(energy == 5500)
+ pars->SetEnergy(AliFMDAnaParameters::k5500);
+
+ if(trigger == 0)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB1);
+ else if(trigger == 1)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB2);
+
+ if(mag==0)
+ pars->SetMagField(AliFMDAnaParameters::k0G);
+ else if(mag==1)
+ pars->SetMagField(AliFMDAnaParameters::k5G);
+
+ if(collsystem == 0)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPP);
+ else if(collsystem == 1)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPbPb);
+
+ pars->PrintStatus();
+
+ std::cout<<"creating background object"<<std::endl;
+ AliFMDAnalysisTaskGenerateCorrection t;
+
+ t.ReadFromFile(filename,store,0);
+ std::cout<<"object created in backgroundFromFile.root "<<std::flush;
+ if(store)
+ std::cout<<" - and stored!"<<std::endl;
+ else
+ std::cout<<" - and not stored!"<<std::endl;
+}
--- /dev/null
+void SubmitFMDCorrections(const Char_t* filename, Bool_t store, Float_t energy, Int_t trigger, Float_t mag, Int_t collsystem, Bool_t realdata=kTRUE) {
+
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libPWG0base");
+ gSystem->Load("libPWG2forward");
+
+ AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
+
+ if(energy == 900)
+ pars->SetEnergy(AliFMDAnaParameters::k900);
+ else if(energy == 7000)
+ pars->SetEnergy(AliFMDAnaParameters::k7000);
+ else if(energy == 10000)
+ pars->SetEnergy(AliFMDAnaParameters::k10000);
+ else if(energy == 14000)
+ pars->SetEnergy(AliFMDAnaParameters::k14000);
+ else if(energy == 5500)
+ pars->SetEnergy(AliFMDAnaParameters::k5500);
+
+ if(trigger == 0)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB1);
+ else if(trigger == 1)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB2);
+
+ if(mag==0)
+ pars->SetMagField(AliFMDAnaParameters::k0G);
+ else if(mag==1)
+ pars->SetMagField(AliFMDAnaParameters::k5G);
+
+ if(collsystem == 0)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPP);
+ else if(collsystem == 1)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPbPb);
+
+ pars->SetRealData(realdata);
+
+ pars->PrintStatus();
+ pars->Init(kTRUE,AliFMDAnaParameters::kBackgroundCorrection);
+ std::cout<<"creating energy distribution object"<<std::endl;
+ AliFMDAnalysisTaskCollector t;
+
+ t.ReadFromFile(filename,store,-1);
+ //std::cout<<"object created in b.root "<<std::flush;
+ if(store)
+ std::cout<<" - and stored!"<<std::endl;
+ else
+ std::cout<<" - and not stored!"<<std::endl;
+}
--- /dev/null
+void SubmitSharingEffCorrection(const Char_t* filename="fmdana.root", Bool_t store, Float_t energy, Int_t trigger, Float_t mag, Int_t collsystem)){
+
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libPWG2forward");
+
+ gStyle->SetTextFont(132);
+ gStyle->SetLabelFont(132,"X");
+ gStyle->SetLabelFont(132,"Y");
+ gStyle->SetLabelFont(132,"Z");
+ gStyle->SetTitleFont(132,"X");
+ gStyle->SetTitleFont(132,"Y");
+ gStyle->SetTitleFont(132,"Z");
+
+ AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
+
+ if(energy == 900)
+ pars->SetEnergy(AliFMDAnaParameters::k900);
+ else if(energy == 7000)
+ pars->SetEnergy(AliFMDAnaParameters::k7000);
+ else if(energy == 10000)
+ pars->SetEnergy(AliFMDAnaParameters::k10000);
+ else if(energy == 14000)
+ pars->SetEnergy(AliFMDAnaParameters::k14000);
+
+ if(trigger == 0)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB1);
+ else if(trigger == 1)
+ pars->SetTriggerDefinition(AliFMDAnaParameters::kMB2);
+
+ if(mag==0)
+ pars->SetMagField(AliFMDAnaParameters::k0G);
+ else if(mag==1)
+ pars->SetMagField(AliFMDAnaParameters::k5G);
+
+ if(collsystem == 0)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPP);
+ else if(collsystem == 1)
+ pars->SetCollisionSystem(AliFMDAnaParameters::kPbPb);
+
+ pars->PrintStatus();
+
+ std::cout<<"creating sharing efficiency object"<<std::endl;
+ AliFMDDndeta t;
+ t.SetNbinsToCut(2);
+ // t.SetVtxCut(2);
+ t.Init(filename);
+
+ t.CreateSharingEfficiency(filename,store);
+ if(store)
+ std::cout<<" - and stored!"<<std::endl;
+ else
+ std::cout<<" - and not stored!"<<std::endl;
+
+}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff