#include "AliESDVertex.h"
#include "TMath.h"
#include "AliFMDAnaParameters.h"
-#include "AliFMDGeometry.h"
+//#include "AliFMDGeometry.h"
#include "AliGenEventHeader.h"
-#include "TDatabasePDG.h"
-#include "TParticlePDG.h"
+#include "AliHeader.h"
+//#include "TDatabasePDG.h"
+//#include "TParticlePDG.h"
+#include "AliFMDStripIndex.h"
ClassImp(AliFMDAnalysisTaskDndeta)
: fDebug(0),
fOutputList(0),
fInputList(0),
- fArray(0),
- fInputArray(0),
fVertexString(0x0),
fNevents(),
fNMCevents(),
fStandalone(kTRUE),
- fMCevent(0)
+ fLastTrackByStrip(0)
{
// Default constructor
DefineInput (0, TList::Class());
fDebug(0),
fOutputList(0),
fInputList(0),
- fArray(),
- fInputArray(0),
fVertexString(0x0),
fNevents(),
fNMCevents(),
fStandalone(kTRUE),
- fMCevent(0)
+ fLastTrackByStrip(0)
{
fStandalone = SE;
if(fStandalone) {
{
AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
- fArray.SetName("FMD");
- fArray.SetOwner();
-
if(!fOutputList)
fOutputList = new TList();
fOutputList->SetName("BackgroundCorrected");
TH2F* hMult = 0;
+ TH1F* hHits = 0;
TH1F* hPrimVertexBin = 0;
- TH2F* hBg = pars->GetBackgroundCorrection(1, 'I', 0);
+ TH2F* hBgTmp = pars->GetBackgroundCorrection(1, 'I', 0);
TH1F* hPrimary = new TH1F("hMultvsEta","hMultvsEta",
- hBg->GetNbinsX(),
- hBg->GetXaxis()->GetXmin(),
- hBg->GetXaxis()->GetXmax());
+ hBgTmp->GetNbinsX(),
+ hBgTmp->GetXaxis()->GetXmin(),
+ hBgTmp->GetXaxis()->GetXmax());
hPrimary->Sumw2();
fOutputList->Add(hPrimary);
Int_t nVtxbins = pars->GetNvtxBins();
for(Int_t det =1; det<=3;det++)
{
- TObjArray* detArray = new TObjArray();
- detArray->SetName(Form("FMD%d",det));
- fArray.AddAtAndExpand(detArray,det);
Int_t nRings = (det==1 ? 1 : 2);
for(Int_t ring = 0;ring<nRings;ring++)
{
Char_t ringChar = (ring == 0 ? 'I' : 'O');
Int_t nSec = (ring == 0 ? 20 : 40);
- TObjArray* vtxArray = new TObjArray();
- vtxArray->SetName(Form("FMD%d%c",det,ringChar));
- detArray->AddAtAndExpand(vtxArray,ring);
+
for(Int_t i = 0; i< nVtxbins; i++) {
TH2F* hBg = pars->GetBackgroundCorrection(det, ringChar, i);
hMult = new TH2F(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i),
hBg->GetXaxis()->GetXmin(),
hBg->GetXaxis()->GetXmax(),
nSec, 0, 2*TMath::Pi());
+
+ hHits = new TH1F(Form("hMCHits_FMD%d%c_vtxbin%d",det,ringChar,i),Form("hMCHits_FMD%d%c_vtxbin%d",det,ringChar,i),
+ hBg->GetNbinsX(),
+ hBg->GetXaxis()->GetXmin(),
+ hBg->GetXaxis()->GetXmax());
+
+
+
hMult->Sumw2();
+ hHits->Sumw2();
fOutputList->Add(hMult);
- vtxArray->AddAtAndExpand(hMult,i);
-
+ fOutputList->Add(hHits);
+
}
}
}
hPrimVertexBin = new TH1F(Form("primmult_vtxbin%d",i),
Form("primmult_vtxbin%d",i),
- hBg->GetNbinsX(),
- hBg->GetXaxis()->GetXmin(),
- hBg->GetXaxis()->GetXmax());
+ hBgTmp->GetNbinsX(),
+ hBgTmp->GetXaxis()->GetXmin(),
+ hBgTmp->GetXaxis()->GetXmax());
hPrimVertexBin->Sumw2();
fOutputList->Add(hPrimVertexBin);
{
if(fStandalone) {
fInputList = (TList*)GetInputData(0);
- fVertexString = (TObjString*)GetInputData(1);
}
}
//_____________________________________________________________________
void AliFMDAnalysisTaskDndeta::Exec(Option_t */*option*/)
{
+ AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
+ fVertexString = (TObjString*)fInputList->At(0);
+
Int_t vtxbin = fVertexString->GetString().Atoi();
+
fNevents.Fill(vtxbin);
+
for(UShort_t det=1;det<=3;det++) {
- //TObjArray* detInputArray = (TObjArray*)fInputArray->At(det);
- TObjArray* detArray = (TObjArray*)fArray.At(det);
Int_t nRings = (det==1 ? 1 : 2);
for (UShort_t ir = 0; ir < nRings; ir++) {
Char_t ringChar = (ir == 0 ? 'I' : 'O');
- //TObjArray* vtxInputArray = (TObjArray*)detInputArray->At(ir);
- TObjArray* vtxArray = (TObjArray*)detArray->At(ir);
- TH2F* hMultTotal = (TH2F*)vtxArray->At(vtxbin);
-
+
+ TH2F* hMultTotal = (TH2F*)fOutputList->FindObject(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
+
TH2F* hMultInput = (TH2F*)fInputList->FindObject(Form("mult_FMD%d%c_vtxbin%d",det,ringChar,vtxbin));
hMultTotal->Add(hMultInput);
-
-
+
}
}
- if(fMCevent)
+ if(pars->GetProcessPrimary())
ProcessPrimary();
if(fStandalone) {
Int_t nVtxbins = pars->GetNvtxBins();
for(UShort_t det=1;det<=3;det++) {
- TObjArray* detArray = (TObjArray*)fArray.At(det);
Int_t nRings = (det==1 ? 1 : 2);
for (UShort_t ir = 0; ir < nRings; ir++) {
- TObjArray* vtxArray = (TObjArray*)detArray->At(ir);
Char_t ringChar = (ir == 0 ? 'I' : 'O');
for(Int_t i =0; i<nVtxbins; i++) {
- TH2F* hMultTotal = (TH2F*)vtxArray->At(i);
+
+ TH2F* hMultTotal = (TH2F*)fOutputList->FindObject(Form("dNdeta_FMD%d%c_vtxbin%d",det,ringChar,i));
+ TH2F* hMultTrVtx = (TH2F*)hMultTotal->Clone(Form("dNdeta_FMD%d%c_TrVtx_vtxbin%d",det,ringChar,i));
+
+ hMultTotal->Scale(pars->GetEventSelectionEfficiency(i));
TH1D* hMultProj = hMultTotal->ProjectionX(Form("dNdeta_FMD%d%c_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY());
+ TH1D* hMultProjTrVtx = hMultTrVtx->ProjectionX(Form("dNdeta_FMD%d%c_TrVtx_vtxbin%d_proj",det,ringChar,i),1,hMultTotal->GetNbinsY());
+ fOutputList->Add(hMultTrVtx);
fOutputList->Add(hMultProj);
+ fOutputList->Add(hMultProjTrVtx);
}
}
}
}
//_____________________________________________________________________
void AliFMDAnalysisTaskDndeta::ProcessPrimary() {
-
+
+ AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
+ AliMCEvent* mcEvent = eventHandler->MCEvent();
+ if(!mcEvent)
+ return;
+
+ fLastTrackByStrip.Reset(-1);
+
AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
AliMCParticle* particle = 0;
- AliStack* stack = fMCevent->Stack();
+ AliStack* stack = mcEvent->Stack();
TH1F* hPrimary = (TH1F*)fOutputList->FindObject("hMultvsEta");
+ AliHeader* header = mcEvent->Header();
+ AliGenEventHeader* genHeader = header->GenEventHeader();
+
+ 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;
Bool_t firstTrack = kTRUE;
- Int_t nTracks = fMCevent->GetNumberOfTracks();
+
+ // 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();
+
for(Int_t i = 0 ;i<nTracks;i++) {
- particle = fMCevent->GetTrack(i);
+ particle = (AliMCParticle*) mcEvent->GetTrack(i);
if(!particle)
continue;
- if(TMath::Abs(particle->Zv()) > pars->GetVtxCutZ())
- continue;
-
+
if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
hPrimary->Fill(particle->Eta());
- Double_t delta = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
- Double_t vertexBinDouble = (particle->Zv() + pars->GetVtxCutZ()) / delta;
- Int_t vertexBin = (Int_t)vertexBinDouble;
+
TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin));
hPrimVtxBin->Fill(particle->Eta());
fNMCevents.Fill(vertexBin);
firstTrack = kFALSE;
}
+
+ }
+ if(pars->GetProcessHits()) {
+ for(Int_t j=0; j<particle->GetNumberOfTrackReferences();j++) {
+
+ AliTrackReference* ref = particle->GetTrackReference(j);
+ UShort_t det,sec,strip;
+ Char_t ring;
+ 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;
+
+ Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex.At(2));//-1*TMath::Log(TMath::Tan(0.5*theta));
+ TH1F* hHits = (TH1F*)fOutputList->FindObject(Form("hMCHits_FMD%d%c_vtxbin%d",det,ring,vertexBin));
+ hHits->Fill(eta);
+ Float_t nstrips = (ring =='O' ? 256 : 512);
+
+ 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;
+
+ }
+ }
}
-
}
-
-
-
-
-
}
//_____________________________________________________________________
//
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff