#include <TH2F.h>
#include <TParticle.h>
#include <TMath.h>
+#include <TParameter.h>
//====================================================================
AliFMDMCEventInspector::AliFMDMCEventInspector()
fHVertex(0),
fHPhiR(0),
fHB(0),
+ fHMcC(0),
fHBvsPart(0),
fHBvsBin(0),
fHBvsCent(0),
fHVzComp(0),
fHCentVsPart(0),
fHCentVsBin(0),
+ fHCentVsMcC(0),
fProduction("")
{
//
fHVertex(0),
fHPhiR(0),
fHB(0),
+ fHMcC(0),
fHBvsPart(0),
fHBvsBin(0),
fHBvsCent(0),
fHVzComp(0),
fHCentVsPart(0),
fHCentVsBin(0),
+ fHCentVsMcC(0),
fProduction("")
{
//
// Parameters:
// name Name of object
//
+ fMC = true;
}
//____________________________________________________________________
fHVertex(0),
fHPhiR(0),
fHB(0),
+ fHMcC(0),
fHBvsPart(0),
fHBvsBin(0),
fHBvsCent(0),
fHVzComp(0),
fHCentVsPart(0),
fHCentVsBin(0),
+ fHCentVsMcC(0),
fProduction("")
{
//
//____________________________________________________________________
void
-AliFMDMCEventInspector::Init(const TAxis& vtxAxis)
+AliFMDMCEventInspector::SetupForData(const TAxis& vtxAxis)
{
//
// Initialize the object
// Parameters:
// vtxAxis Vertex axis in use
//
- AliFMDEventInspector::Init(vtxAxis);
+
+ // We temporary disable the displaced vertices so we can initialize
+ // the routine ourselves.
+ Bool_t saveDisplaced = AllowDisplaced();
+ if (saveDisplaced) SetVertexMethod(kNormal);
+ AliFMDEventInspector::SetupForData(vtxAxis);
+ if (saveDisplaced) SetVertexMethod(kDisplaced);
Int_t maxPart = 450;
Int_t maxBin = 225;
fHB->SetDirectory(0);
fList->Add(fHB);
+ fHMcC = static_cast<TH1F*>(fHCent->Clone("mcC"));
+ fHMcC->SetFillColor(kCyan+2);
+ fHMcC->SetDirectory(0);
+ fList->Add(fHMcC);
+
fHBvsPart = new TH2F("bVsParticipants", "Impact parameter vs Participants",
5*maxB, 0, maxB, maxPart, -.5, maxPart-.5);
fHBvsPart->SetXTitle("b [fm]");
fHCentVsBin->SetZTitle("Event");
fHCentVsBin->SetDirectory(0);
fList->Add(fHCentVsBin);
-}
-//____________________________________________________________________
-void
-AliFMDMCEventInspector::StoreInformation(Int_t runNo)
-{
- // Store information about running conditions in the output list
- if (!fList) return;
- AliFMDEventInspector::StoreInformation(runNo);
- TNamed* mc = new TNamed("mc", fProduction.Data());
- mc->SetUniqueID(1);
- fList->Add(mc);
+ Int_t nC = fHCent->GetNbinsX();
+ Double_t cL = fHCent->GetXaxis()->GetXmin();
+ Double_t cH = fHCent->GetXaxis()->GetXmax();
+ fHCentVsMcC = new TH2F("centralityRecoVsMC",
+ "Centrality from reconstruction vs MC derived",
+ nC, cL, cH, nC, cL, cH);
+ fHCentVsMcC->SetDirectory(0);
+ fHCentVsMcC->SetStats(0);
+ fHCentVsMcC->SetXTitle("Centralty from Reco [%]");
+ fHCentVsMcC->SetYTitle("Centralty derived from Impact Par. [%]");
+ fHCentVsMcC->SetZTitle("Events");
+ fList->Add(fHCentVsMcC);
+
+ if (AllowDisplaced()) fDisplacedVertex.SetupForData(fList, "", true);
}
namespace
UShort_t& ivz,
Double_t& vz,
Double_t& b,
+ Double_t& c,
Int_t& npart,
Int_t& nbin,
Double_t& phiR)
// means outside of the defined vertex range
// vz On return, the z position of the interaction
// b On return, the impact parameter - < 0 if not available
+ // c On return, centrality estimate - < 0 if not available
// phiR On return, the event plane angle - < 0 if not available
//
// Return:
npart = 0;
nbin = 0;
b = -1;
+ c = -1;
+ phiR = -1111;
if (colGeometry) {
b = colGeometry->ImpactParameter();
phi = colGeometry->ReactionPlaneAngle();
colGeometry->TargetParticipants());
nbin = colGeometry->NN();
}
+ if (fDebug && !colGeometry) {
+ AliWarningF("Collision header of class %s is not a CollisionHeader",
+ genHeader->ClassName());
+ }
+
if(pythiaHeader) {
Int_t pythiaType = pythiaHeader->ProcessType();
// 92 and 93 are SD
npart = 2; // Always 2 protons
nbin = 1; // Always 1 binary collision
}
+ if (b >= 0) {
+ if (b < 3.5) c = 2.5; //0-5%
+ else if (b < 4.95) c = 7.5; //5-10%
+ else if (b < 6.98) c = 15; //10-20%
+ else if (b < 8.55) c = 25; //20-30%
+ else if (b < 9.88) c = 35; //30-40%
+ else if (b < 11.04) c = 45; //40-50%
+ else c = 55; //50-60%
+ }
if(dpmHeader) { // Also an AliCollisionGeometry
Int_t processType = dpmHeader->ProcessType();
// 1 & 4 are ND
else if (phi > 2*TMath::Pi()) phi -= 2*TMath::Pi();
else break;
}
+ phiR = phi;
}
// Do a check on particles
genHeader->PrimaryVertex(vtx);
vz = vtx[2];
+ DMSG(fDebug, 2, "vz=%f, phiR=%f, b=%f, npart=%d, nbin=%d",
+ vz, phiR, b, npart, nbin);
+
fHVertex->Fill(vz);
fHPhiR->Fill(phiR);
fHB->Fill(b);
+ fHMcC->Fill(c);
fHBvsPart->Fill(b, npart);
fHBvsBin->Fill(b, nbin);
+ if(AllowDisplaced()) {
+#if 0
+ // Put the vertex at fixed locations
+ Double_t zvtx = vz;
+ Double_t ratio = zvtx/37.5;
+ if(ratio > 0) ratio = ratio + 0.5;
+ if(ratio < 0) ratio = ratio - 0.5;
+ Int_t ratioInt = Int_t(ratio);
+ zvtx = 37.5*((Double_t)ratioInt);
+ if(TMath::Abs(zvtx) > 999)
+ return kBadVertex;
+#endif
+ if (!fDisplacedVertex.ProcessMC(event))
+ return kBadVertex;
+ if (fDisplacedVertex.IsSatellite())
+ vz = fDisplacedVertex.GetVertexZ();
+ }
+
// Check for the vertex bin
ivz = fVtxAxis.FindBin(vz);
if (ivz <= 0 || ivz > fHEventsTr->GetXaxis()->GetNbins()) {
return kOk;
}
+//____________________________________________________________________
+Bool_t
+AliFMDMCEventInspector::ReadCentrality(const AliESDEvent& esd,
+ Double_t& cent,
+ UShort_t& qual) const
+{
+ //
+ // Read centrality from event
+ //
+ // Parameters:
+ // esd Event
+ // cent On return, the centrality or negative if not found
+ //
+ // Return:
+ // False on error, true otherwise
+ //
+ Bool_t ret = AliFMDEventInspector::ReadCentrality(esd, cent, qual);
+ if (qual != 0) {
+ AliCentrality* centObj = const_cast<AliESDEvent&>(esd).GetCentrality();
+ if (!centObj) return ret;
+
+ // For MC, we allow `bad' centrality selections
+ cent = centObj->GetCentralityPercentileUnchecked(fCentMethod);
+ }
+ return ret;
+}
//____________________________________________________________________
namespace {
return false;
// Rapidity shift
- Double_t m02s = 1 - 2 * p1->Energy() / fEnergy;
- Double_t m12s = 1 - 2 * p2->Energy() / fEnergy;
+ Double_t m02s = (fEnergy > 0 ? 1 - 2 * p1->Energy() / fEnergy : 0);
+ Double_t m12s = (fEnergy > 0 ? 1 - 2 * p2->Energy() / fEnergy : 0);
if (arm == 0 && m02s > xiMin && m02s < xiMax) return true;
if (arm == 1 && m12s > xiMin && m12s < xiMax) return true;
return false;
}
-//____________________________________________________________________
-Bool_t
-AliFMDMCEventInspector::ReadCentrality(const AliESDEvent* esd,
- Double_t& cent,
- UShort_t& qual) const
-{
- //
- // Read centrality from event
- //
- // Parameters:
- // esd Event
- // cent On return, the centrality or negative if not found
- //
- // Return:
- // False on error, true otherwise
- //
- cent = -1;
- qual = 0;
- AliCentrality* centObj = const_cast<AliESDEvent*>(esd)->GetCentrality();
- if (!centObj) return true;
-AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
- Bool_t isMC = am->GetMCtruthEventHandler() != 0;
-
- //std::cout<<fUseDisplacedVertices<<" "<<isMC<<std::endl;
-
- if(fUseDisplacedVertices && isMC) {
-
-
- AliMCEventHandler* mchandler = static_cast<AliMCEventHandler*>(am->GetMCtruthEventHandler());
- AliMCEvent* mcevent = mchandler->MCEvent();
-
- AliHeader* header = mcevent->Header();
- AliGenEventHeader* genHeader = header->GenEventHeader();
- AliCollisionGeometry* colGeometry =
- dynamic_cast<AliCollisionGeometry*>(genHeader);
- Double_t b = -1;
- if (colGeometry)
- b = colGeometry->ImpactParameter();
- //std::cout<<"Hallo!! "<<b<<std::endl;
- cent = -1;
- if(b<3.5 && b >0) cent = 2.5; //0-5%
- if(b>3.5 && b<4.95) cent = 7.5; //5-10%
- if(b>4.95 && b<6.98) cent = 15; //10-20%
- if(b>6.98 && b<8.55) cent = 25; //20-30%
- if(b>8.55 && b<9.88) cent = 35; //30-40%
- if(b>9.88 && b<11.04) cent = 45; //40-50%
- if(b>11.04) cent = 55; //50-60%
- //cent = 10;
- qual = 0;
- }
- else {
-
- qual = centObj->GetQuality();
- if (qual == 0x8) // Ignore ZDC outliers
- cent = centObj->GetCentralityPercentileUnchecked("V0M");
- else
- cent = centObj->GetCentralityPercentile("V0M");
- }
- return true;
-}
//____________________________________________________________________
Bool_t
AliFMDMCEventInspector::CompareResults(Double_t vz, Double_t trueVz,
- Double_t cent, Double_t b,
+ Double_t cent, Double_t mcC,
+ Double_t b,
Int_t npart, Int_t nbin)
{
fHVzComp->Fill(trueVz, vz);
fHBvsCent->Fill(b, cent);
fHCentVsPart->Fill(npart, cent);
fHCentVsBin->Fill(nbin, cent);
+ fHCentVsMcC->Fill(cent, mcC);
return true;
}
+
+
//
// EOF
//
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff