#include "AliFMDStripIndex.h"
#include "AliMCEvent.h"
// #include "AliFMDAnaParameters.h"
+#include "AliESDFMD.h"
#include "AliLog.h"
#include <TH2D.h>
+#include <TProfile2D.h>
ClassImp(AliFMDMCDensityCalculator)
#if 0
; // For Emacs
#endif
+//____________________________________________________________________
+AliFMDMCDensityCalculator::~AliFMDMCDensityCalculator()
+{
+ if (fComps) fComps->Clear();
+ if (fFMD1i) delete fFMD1i;
+ if (fFMD2i) delete fFMD2i;
+ if (fFMD2o) delete fFMD2o;
+ if (fFMD3i) delete fFMD3i;
+ if (fFMD3o) delete fFMD3o;
+ if (fFMD1iC) delete fFMD1iC;
+ if (fFMD2iC) delete fFMD2iC;
+ if (fFMD2oC) delete fFMD2oC;
+ if (fFMD3iC) delete fFMD3iC;
+ if (fFMD3oC) delete fFMD3oC;
+}
//____________________________________________________________________
AliFMDMCDensityCalculator&
}
+
//____________________________________________________________________
Bool_t
-AliFMDMCDensityCalculator::CalculateMC(const AliMCEvent& event,
- AliForwardUtil::Histos& hists,
- Double_t vz,
- UShort_t vtxbin)
+AliFMDMCDensityCalculator::CalculateMC(const AliESDFMD& fmd,
+ AliForwardUtil::Histos& hists)
{
- Int_t nTracks = event.GetNumberOfTracks();
- for (Int_t iTr = 0; iTr < nTracks; iTr++) {
- AliMCParticle* particle =
- static_cast<AliMCParticle*>(event.GetTrack(iTr));
-
- // Check the returned particle
- if (!particle) continue;
-
- // Check if this charged and a primary
- Bool_t isCharged = particle->Charge() != 0;
- if (!isCharged) continue;
+ for (UShort_t d=1; d<=3; d++) {
+ UShort_t nr = (d == 1 ? 1 : 2);
+ for (UShort_t q=0; q<nr; q++) {
+ Char_t r = (q == 0 ? 'I' : 'O');
+ UShort_t ns= (q == 0 ? 20 : 40);
+ UShort_t nt= (q == 0 ? 512 : 256);
+ TH2D* h = hists.Get(d,r);
- Int_t nTrRef = particle->GetNumberOfTrackReferences();
- for (Int_t iTrRef = 0; iTrRef < nTrRef; iTrRef++) {
- AliTrackReference* ref = particle->GetTrackReference(iTrRef);
-
- // Check existence
- if (!ref) continue;
-
- // Check that we hit an FMD element
- if (ref->DetectorId() != AliTrackReference::kFMD)
- continue;
-
- // Get the detector coordinates
- UShort_t d, s, t;
- Char_t r;
- AliFMDStripIndex::Unpack(ref->UserId(), d, r, s, t);
-
- Double_t x = ref->X();
- Double_t y = ref->Y();
- Double_t z = ref->Z()-vz;
- Double_t rr = TMath::Sqrt(x*x+y*y);
- Double_t phi = TMath::ATan2(y,x);
- Double_t theta= TMath::ATan2(rr,z);
- Double_t eta = -TMath::Log(TMath::Tan(theta/2));
-
- Float_t c = Correction(d,r,s,t,vtxbin,eta,false);
- fCorrections->Fill(c);
-
- TH2D* h = hists.Get(d,r);
- h->Fill(eta,phi, 1 * c);
+ for (UShort_t s=0; s<ns; s++) {
+ for (UShort_t t=0; t<nt; t++) {
+ Float_t mult = fmd.Multiplicity(d,r,s,t);
+
+ if (mult == 0 || mult > 20) continue;
+
+ Float_t phi = fmd.Phi(d,r,s,t) / 180 * TMath::Pi();
+ Float_t eta = fmd.Eta(d,r,s,t);
+
+ Float_t corr = Correction(d, r, s, t, eta, 0, false);
+ Float_t sig = (corr <= 0 ? 0 : mult / corr);
+ h->Fill(eta,phi,sig);
+ }
+ }
}
}
return kTRUE;
}
+//____________________________________________________________________
+void
+AliFMDMCDensityCalculator::Init(const TAxis& eAxis)
+{
+ fFMD1i = Make(1,'I',eAxis);
+ fFMD2i = Make(2,'I',eAxis);
+ fFMD2o = Make(2,'O',eAxis);
+ fFMD3i = Make(3,'I',eAxis);
+ fFMD3o = Make(3,'O',eAxis);
+ fFMD1iC = Make(1,'I');
+ fFMD2iC = Make(2,'I');
+ fFMD2oC = Make(2,'O');
+ fFMD3iC = Make(3,'I');
+ fFMD3oC = Make(3,'O');
+
+ fComps->Add(fFMD1i);
+ fComps->Add(fFMD2i);
+ fComps->Add(fFMD2o);
+ fComps->Add(fFMD3i);
+ fComps->Add(fFMD3o);
+ fComps->Add(fFMD1iC);
+ fComps->Add(fFMD2iC);
+ fComps->Add(fFMD2oC);
+ fComps->Add(fFMD3iC);
+ fComps->Add(fFMD3oC);
+}
+
+//____________________________________________________________________
+TProfile2D*
+AliFMDMCDensityCalculator::Make(UShort_t d, Char_t r,
+ const TAxis& axis) const
+{
+ TProfile2D* ret = new TProfile2D(Form("FMD%d%c_esd_vs_mc", d, r),
+ Form("ESD/MC signal for FMD%d%c", d, r),
+ axis.GetNbins(),
+ axis.GetXmin(),
+ axis.GetXmax(),
+ (r == 'I' || r == 'i') ? 20 : 40,
+ 0, 2*TMath::Pi());
+ ret->GetXaxis()->SetTitle("#eta");
+ ret->GetYaxis()->SetTitle("#varphi [degrees]");
+ ret->GetZaxis()->SetTitle("#LT incusive density ESD/MC#GT");
+ ret->SetDirectory(0);
+ return ret;
+}
+//____________________________________________________________________
+TH2D*
+AliFMDMCDensityCalculator::Make(UShort_t d, Char_t r) const
+{
+ TH2D* ret = new TH2D(Form("FMD%d%c_corr_mc_esd", d, r),
+ Form("ESD-MC correlation for FMD%d%c", d, r),
+ 200, 0, 20, 200, 0, 20);
+ ret->GetXaxis()->SetTitle("m_{incl} (MC)");
+ ret->GetYaxis()->SetTitle("#Delta/#Delta_{mp} (ESD)");
+ ret->GetZaxis()->SetTitle("Correlation");
+ ret->SetDirectory(0);
+ return ret;
+}
+//____________________________________________________________________
+void
+AliFMDMCDensityCalculator::Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc)
+{
+ if (!esd || !mc) return;
+ TProfile2D* p = 0;
+ TH2D* h = 0;
+ switch (d) {
+ case 1:
+ p = fFMD1i;
+ h = fFMD1iC;
+ break;
+ case 2:
+ p = (r == 'I' || r == 'i' ? fFMD2i : fFMD2o);
+ h = (r == 'I' || r == 'i' ? fFMD2iC : fFMD2oC);
+ break;
+ case 3:
+ p = (r == 'I' || r == 'i' ? fFMD3i : fFMD3o);
+ h = (r == 'I' || r == 'i' ? fFMD3iC : fFMD3oC);
+ break;
+ }
+ if (!p) return;
+
+ for (Int_t iEta = 1; iEta <= esd->GetNbinsX(); iEta++) {
+ Double_t eta = esd->GetXaxis()->GetBinCenter(iEta);
+ for (Int_t iPhi = 1; iPhi <= esd->GetNbinsY(); iPhi++) {
+ Double_t phi = esd->GetYaxis()->GetBinCenter(iPhi);
+ Double_t mEsd = esd->GetBinContent(iEta,iPhi);
+ Double_t mMc = mc->GetBinContent(iEta,iPhi);
+
+ p->Fill(eta, phi, (mMc > 0 ? mEsd / mMc : 0));
+ h->Fill(mMc,mEsd);
+ }
+ }
+}
+
//____________________________________________________________________
Bool_t
-AliFMDMCDensityCalculator::Calculate(const AliESDFMD&,
- AliForwardUtil::Histos&,
- UShort_t,
- Bool_t)
+AliFMDMCDensityCalculator::CompareResults(AliForwardUtil::Histos& esd,
+ AliForwardUtil::Histos& mc)
{
- AliWarning("Method Calculate disabled for this class. If you need this, "
- "make an AliFMDDensityCalculator object instead");
- return kFALSE;
+ Fill(1, 'I', esd.Get(1,'I'), mc.Get(1,'I'));
+ Fill(2, 'I', esd.Get(2,'I'), mc.Get(2,'I'));
+ Fill(2, 'O', esd.Get(2,'O'), mc.Get(2,'O'));
+ Fill(3, 'I', esd.Get(3,'I'), mc.Get(3,'I'));
+ Fill(3, 'O', esd.Get(3,'O'), mc.Get(3,'O'));
+
+ return kTRUE;
}
+//____________________________________________________________________
+void
+AliFMDMCDensityCalculator::DefineOutput(TList* dir)
+{
+ AliFMDDensityCalculator::DefineOutput(dir);
+ TList* d = static_cast<TList*>(dir->FindObject(GetName()));
+
+ fComps = new TList;
+ fComps->SetName("esd_mc_comparison");
+ d->Add(fComps);
+
+}
//____________________________________________________________________
//
// EOF
-#ifndef ALIROOT_PWG2_FORWARD_ANALYSIS2_ALIFMDMCDENSITYCALCULATOR_H
-#define ALIROOT_PWG2_FORWARD_ANALYSIS2_ALIFMDMCDENSITYCALCULATOR_H
+#ifndef ALIFMDMCDENSITYCALCULATOR_H
+#define ALIFMDMCDENSITYCALCULATOR_H
#include "AliFMDDensityCalculator.h"
#include <TList.h>
#include "AliForwardUtil.h"
class AliMCEvent;
+class TH2;
class TH2D;
-class TH1D;
+class TProfile2D;
/**
* This class calculates the inclusive charged particle density
* @par Corrections used:
* - None
*
- * @ingroup pwg2_forward
+ * @ingroup pwg2_forward_algo
+ * @ingroup pwg2_forward_mc
*/
class AliFMDMCDensityCalculator : public AliFMDDensityCalculator
{
/**
* Constructor
*/
- AliFMDMCDensityCalculator() : AliFMDDensityCalculator() {}
+ AliFMDMCDensityCalculator()
+ : AliFMDDensityCalculator(),
+ fFMD1i(0),
+ fFMD2i(0),
+ fFMD2o(0),
+ fFMD3i(0),
+ fFMD3o(0),
+ fFMD1iC(0),
+ fFMD2iC(0),
+ fFMD2oC(0),
+ fFMD3iC(0),
+ fFMD3oC(0),
+ fComps(0)
+ {}
/**
* Constructor
*
* @param name Name of object
*/
AliFMDMCDensityCalculator(const char* name)
- : AliFMDDensityCalculator(name)
+ : AliFMDDensityCalculator(name),
+ fFMD1i(0),
+ fFMD2i(0),
+ fFMD2o(0),
+ fFMD3i(0),
+ fFMD3o(0),
+ fFMD1iC(0),
+ fFMD2iC(0),
+ fFMD2oC(0),
+ fFMD3iC(0),
+ fFMD3oC(0),
+ fComps(0)
{}
/**
* Copy constructor
* @param o Object to copy from
*/
AliFMDMCDensityCalculator(const AliFMDMCDensityCalculator& o)
- : AliFMDDensityCalculator(o)
+ : AliFMDDensityCalculator(o) ,
+ fFMD1i(o.fFMD1i),
+ fFMD2i(o.fFMD2i),
+ fFMD2o(o.fFMD2o),
+ fFMD3i(o.fFMD3i),
+ fFMD3o(o.fFMD3o),
+ fFMD1iC(o.fFMD1iC),
+ fFMD2iC(o.fFMD2iC),
+ fFMD2oC(o.fFMD2oC),
+ fFMD3iC(o.fFMD3iC),
+ fFMD3oC(o.fFMD3oC),
+ fComps(0)
{}
/**
* Destructor
*/
- virtual ~AliFMDMCDensityCalculator() {}
+ virtual ~AliFMDMCDensityCalculator();
/**
* Assignement operator
*
*/
AliFMDMCDensityCalculator& operator=(const AliFMDMCDensityCalculator& o);
/**
- * Do the calculations
- *
- * @param fmd AliESDFMD object (possibly) corrected for sharing
- * @param hists Histogram cache
- * @param vtxBin Vertex bin
- * @param lowFlux Low flux flag.
+ * Initialize this object
*
- * @return true on successs
+ * @param etaAxis Eta axis to use
*/
- virtual Bool_t Calculate(const AliESDFMD& fmd,
- AliForwardUtil::Histos& hists,
- UShort_t vtxBin, Bool_t lowFlux);
+ void Init(const TAxis& etaAxis);
/**
* Calculate the charged particle density from the MC track references.
*
*
* @return true on success
*/
- virtual Bool_t CalculateMC(const AliMCEvent& event,
- AliForwardUtil::Histos& hists,
- Double_t vz,
- UShort_t vtxBin);
+ virtual Bool_t CalculateMC(const AliESDFMD& fmd,
+ AliForwardUtil::Histos& hists);
+
+ /**
+ * Compare the result of analysing the ESD for
+ * the inclusive charged particle density to analysing
+ * MC truth
+ *
+ * @param esd
+ * @param mc
+ *
+ * @return
+ */
+ virtual Bool_t CompareResults(AliForwardUtil::Histos& esd,
+ AliForwardUtil::Histos& mc);
+ /**
+ * Output diagnostic histograms to directory
+ *
+ * @param dir List to write in
+ */
+ void DefineOutput(TList* dir);
protected:
+ /**
+ * MAke comparison profiles
+ *
+ * @param d Detector
+ * @param r Ring
+ * @param axis Eta axis
+ *
+ * @return Newly allocated profile object
+ */
+ TProfile2D* Make(UShort_t d, Char_t r, const TAxis& axis) const;
+ /**
+ * MAke comparison profiles
+ *
+ * @param d Detector
+ * @param r Ring
+ *
+ * @return Newly allocated profile object
+ */
+ TH2D* Make(UShort_t d, Char_t r) const;
+ /**
+ * Fill comparison profiles
+ *
+ * @param d Detector
+ * @param r Ring
+ * @param esd ESD histogram
+ * @param mc MC histogram
+ */
+ void Fill(UShort_t d, Char_t r, TH2* esd, TH2* mc);
+
+ TProfile2D* fFMD1i; // Comparison
+ TProfile2D* fFMD2i; // Comparison
+ TProfile2D* fFMD2o; // Comparison
+ TProfile2D* fFMD3i; // Comparison
+ TProfile2D* fFMD3o; // Comparison
+ TH2D* fFMD1iC; // Correlation in FMD1i
+ TH2D* fFMD2iC; // Correlation in FMD2i
+ TH2D* fFMD2oC; // Correlation in FMD2o
+ TH2D* fFMD3iC; // Correlation in FMD3i
+ TH2D* fFMD3oC; // Correlation in FMD3o
+ TList* fComps; // List of comparisons
ClassDef(AliFMDMCDensityCalculator,1); // Calculate Nch density
};