AliFMDEventInspector, AliFMDDensityCalculator, AliFMDSharingFilter:
authorcholm <cholm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 13 May 2011 14:50:50 +0000 (14:50 +0000)
committercholm <cholm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 13 May 2011 14:50:50 +0000 (14:50 +0000)
  Store information on settings in output file
analysis2/TrainSetup.C: Fixes
AliAODForwardMult AliFMDEventInspector - revert to standard use of
  physics selection - no more AC and E triggers.
AliFMDCorrELossFit utility to get plots for comparions
OtherData - use 2360 for 2760

19 files changed:
PWG2/FORWARD/analysis2/AliAODForwardMult.cxx
PWG2/FORWARD/analysis2/AliFMDCorrELossFit.cxx
PWG2/FORWARD/analysis2/AliFMDCorrELossFit.h
PWG2/FORWARD/analysis2/AliFMDDensityCalculator.cxx
PWG2/FORWARD/analysis2/AliFMDDensityCalculator.h
PWG2/FORWARD/analysis2/AliFMDEnergyFitter.cxx
PWG2/FORWARD/analysis2/AliFMDEnergyFitterTask.cxx
PWG2/FORWARD/analysis2/AliFMDEventInspector.cxx
PWG2/FORWARD/analysis2/AliFMDEventInspector.h
PWG2/FORWARD/analysis2/AliFMDMCEventInspector.cxx
PWG2/FORWARD/analysis2/AliFMDMCEventInspector.h
PWG2/FORWARD/analysis2/AliFMDSharingFilter.cxx
PWG2/FORWARD/analysis2/AliForwardMCMultiplicityTask.cxx
PWG2/FORWARD/analysis2/AliForwardMultiplicityTask.cxx
PWG2/FORWARD/analysis2/AliMCTruthdNdetaTask.cxx
PWG2/FORWARD/analysis2/ForwardAODConfig.C
PWG2/FORWARD/analysis2/MakeAOD.C
PWG2/FORWARD/analysis2/OtherData.C
PWG2/FORWARD/analysis2/TrainSetup.C

index fafa2fb..65505b0 100644 (file)
@@ -316,7 +316,8 @@ AliAODForwardMult::CheckEvent(Int_t    triggerMask,
     if (IsTriggerBits(kNClusterGt0))    hist->AddBinContent(kBinNClusterGt0);
   }
   // Check if we have an event of interest. 
-  if (!IsTriggerBits(triggerMask|kB)) return false;
+  Int_t mask = triggerMask; //|kB
+  if (!IsTriggerBits(mask)) return false;
   
   // Check for pileup
   if (IsTriggerBits(kPileUp)) return false;
index 3991420..76117b4 100644 (file)
@@ -964,39 +964,12 @@ namespace {
 #define IDX2RING(I) (i == 0 || i == 1 || i == 3 ? 'I' : 'O')
 #define IDX2DET(I)  (i == 0 ? 1 : (i == 1 || i == 2 ? 2 : 3))
 //____________________________________________________________________
-void
-AliFMDCorrELossFit::Draw(Option_t* option)
+TList*
+AliFMDCorrELossFit::GetStacks(Bool_t err, Bool_t rel, UShort_t maxN) const
 {
-  // 
-  // Draw this object 
-  // 
-  // Parameters:
-  //    option Options.  Possible values are 
-  //  - err Plot error bars 
-  //
-  TString opt(Form("nostack %s", option));
-  opt.ToLower();
-  Bool_t  rel = (opt.Contains("relative"));
-  Bool_t  err = (opt.Contains("error"));
-  if (rel) opt.ReplaceAll("relative","");
-  if (err) opt.ReplaceAll("error","");
-  Int_t    nRings = fRings.GetEntriesFast();
-  UShort_t maxN   = 0;
-  for (Int_t i = 0; i < nRings; i++) { 
-    if (!fRings.At(i)) continue;
-    TObjArray* a     = static_cast<TObjArray*>(fRings.At(i));
-    Int_t      nFits = a->GetEntriesFast();
-
-    for (Int_t j = 0; j < nFits; j++) {
-      ELossFit* fit = static_cast<ELossFit*>(a->At(j));
-      if (!fit) continue;
-      maxN          = TMath::Max(maxN, UShort_t(fit->fN));
-    }
-  }
-  // AliInfo(Form("Maximum N is %d", maxN));
-  Int_t nPad = 6+maxN-1; // 7 regular params, and maxN-1 weights
-  TVirtualPad* pad = gPad;
-  pad->Divide(2, (nPad+1)/2, 0.1, 0, 0);
+  // Get a list of THStacks 
+  Int_t nRings = fRings.GetEntriesFast();
+  Int_t nPad   = 6+maxN-1; // 7 regular params, and maxN-1 weights
 
   enum { 
     kChi2nu = 0, 
@@ -1014,18 +987,18 @@ AliFMDCorrELossFit::Draw(Option_t* option)
   THStack* sSigma;
   // THStack* sigman;
   THStack* n;
-  TList stacks;
-  stacks.AddAt(sChi2nu= new THStack("chi2",   "#chi^{2}/#nu"), kChi2nu);
-  stacks.AddAt(sC     = new THStack("c",       "C"),           kC);
-  stacks.AddAt(sDelta = new THStack("delta",  "#Delta_{mp}"),  kDelta);
-  stacks.AddAt(sXi    = new THStack("xi",     "#xi"),          kXi);
-  stacks.AddAt(sSigma = new THStack("sigma",  "#sigma"),       kSigma);
-  //stacks.AddAt(sigman= new THStack("sigman", "#sigma_{n}"),   5);
-  stacks.AddAt(n     = new THStack("n",      "N"),            kN);
+  TList* stacks = new TList;
+  stacks->AddAt(sChi2nu= new THStack("chi2",   "#chi^{2}/#nu"), kChi2nu);
+  stacks->AddAt(sC     = new THStack("c",       "C"),           kC);
+  stacks->AddAt(sDelta = new THStack("delta",  "#Delta_{mp}"),  kDelta);
+  stacks->AddAt(sXi    = new THStack("xi",     "#xi"),          kXi);
+  stacks->AddAt(sSigma = new THStack("sigma",  "#sigma"),       kSigma);
+  //stacks->AddAt(sigman= new THStack("sigman", "#sigma_{n}"),   5);
+  stacks->AddAt(n     = new THStack("n",      "N"),            kN);
   for (Int_t i = 1; i <= maxN; i++) {
-    stacks.AddAt(new THStack(Form("a_%02d", i+1), Form("a_{%d}", i+1)), kN+i);
+    stacks->AddAt(new THStack(Form("a_%02d", i+1), Form("a_{%d}", i+1)), kN+i);
   }
-
+  
   TArrayD min(nPad);
   TArrayD max(nPad);
   min.Reset(100000);
@@ -1060,7 +1033,7 @@ AliFMDCorrELossFit::Draw(Option_t* option)
 
     for (Int_t k = 1; k <= maxN; k++) { 
       hA[k-1] = MakeHist(fEtaAxis,a->GetName(), Form("a%02d",k+1), color);
-      static_cast<THStack*>(stacks.At(kN+k))->Add(hA[k-1], ho);
+      static_cast<THStack*>(stacks->At(kN+k))->Add(hA[k-1], ho);
     }
                          
     for (Int_t j = 0; j < nFits; j++) {
@@ -1130,6 +1103,47 @@ AliFMDCorrELossFit::Draw(Option_t* option)
       }
     }
   }
+  return stacks;
+}
+
+//____________________________________________________________________
+void
+AliFMDCorrELossFit::Draw(Option_t* option)
+{
+  // 
+  // Draw this object 
+  // 
+  // Parameters:
+  //    option Options.  Possible values are 
+  //  - err Plot error bars 
+  //
+  TString opt(Form("nostack %s", option));
+  opt.ToLower();
+  Bool_t  rel = (opt.Contains("relative"));
+  Bool_t  err = (opt.Contains("error"));
+  if (rel) opt.ReplaceAll("relative","");
+  if (err) opt.ReplaceAll("error","");
+
+  UShort_t maxN   = 0;
+  Int_t nRings = fRings.GetEntriesFast();
+  for (Int_t i = 0; i < nRings; i++) { 
+    if (!fRings.At(i)) continue;
+    TObjArray* a     = static_cast<TObjArray*>(fRings.At(i));
+    Int_t      nFits = a->GetEntriesFast();
+
+    for (Int_t j = 0; j < nFits; j++) {
+      ELossFit* fit = static_cast<ELossFit*>(a->At(j));
+      if (!fit) continue;
+      maxN          = TMath::Max(maxN, UShort_t(fit->fN));
+    }
+  }
+  // AliInfo(Form("Maximum N is %d", maxN));
+  Int_t nPad = 6+maxN-1; // 7 regular params, and maxN-1 weights
+  TVirtualPad* pad = gPad;
+  pad->Divide(2, (nPad+1)/2, 0.1, 0, 0);
+
+  TList* stacks = GetStacks(err, rel, maxN);
+
   Int_t nPad2 = (nPad+1) / 2;
   for (Int_t i = 0; i < nPad; i++) {
     Int_t iPad = 1 + i/nPad2 + 2 * (i % nPad2);
@@ -1141,11 +1155,13 @@ AliFMDCorrELossFit::Draw(Option_t* option)
     p->SetGridy();
     if (rel && i != 0 && i != 6 && i != 5 && i != 4) p->SetLogy();
 
-    Double_t powMax = TMath::Log10(max[i]);
-    Double_t powMin = min[i] <= 0 ? powMax : TMath::Log10(min[i]);
-    if (powMax-powMin > 2. && min[i] != 0) p->SetLogy();
 
-    THStack* stack = static_cast<THStack*>(stacks.At(i));
+    THStack* stack = static_cast<THStack*>(stacks->At(i));
+
+    // Double_t powMax = TMath::Log10(max[i]);
+    // Double_t powMin = min[i] <= 0 ? powMax : TMath::Log10(min[i]);
+    // if (powMax-powMin > 2. && min[i] != 0) p->SetLogy();
+
     // stack->SetMinimum(min[i]);
     // stack->SetMaximum(max[i]);
     stack->Draw(opt.Data());
index 185ad3f..8eb83c5 100644 (file)
@@ -612,6 +612,15 @@ public:
    *
    */
   void Print(Option_t* option="R") const; //*MENU*
+    /** 
+     * Get a list of THStack - one for each parameter 
+     * 
+     * @param err Show errors
+     * @param rel Show relative errors 
+     * 
+     * @return List of THStack
+     */
+  TList* GetStacks(Bool_t err, Bool_t rel, UShort_t maxN=5) const;
   /* @} */
 protected:
   /** 
index 283cb61..ddbabcd 100644 (file)
@@ -14,6 +14,7 @@
 #include <TProfile.h>
 #include <THStack.h>
 #include <TROOT.h>
+#include <TParameter.h>
 #include <iostream>
 #include <iomanip>
 
@@ -118,7 +119,6 @@ AliFMDDensityCalculator::AliFMDDensityCalculator(const char* title)
   fLowCuts->SetXTitle("#eta");
   fLowCuts->SetDirectory(0);
 
-  for (Int_t i = 0; i < 5; i++) fMultCuts[i] = 0;
 }
 
 //____________________________________________________________________
@@ -261,11 +261,11 @@ AliFMDDensityCalculator::SetMultCuts(Double_t fmd1i,
                                     Double_t fmd3i, 
                                     Double_t fmd3o) 
 {
-  fMultCuts[0] = fmd1i;
-  fMultCuts[1] = fmd2i;
-  fMultCuts[2] = fmd2o;
-  fMultCuts[3] = fmd3i;
-  fMultCuts[4] = fmd3o;
+  GetRingHistos(1,'I')->fMultCut = fmd1i;
+  GetRingHistos(2,'I')->fMultCut = fmd2i;
+  GetRingHistos(2,'O')->fMultCut = fmd2o;
+  GetRingHistos(3,'I')->fMultCut = fmd3i;
+  GetRingHistos(3,'O')->fMultCut = fmd3o;
   fMultCut = (fmd1i+fmd2i+fmd2o+fmd3i+fmd3o) / 5;
 }
 
@@ -282,9 +282,10 @@ AliFMDDensityCalculator::GetMultCut(UShort_t d, Char_t r, Int_t ieta,
   // Return:
   //    Lower cut on multiplicity
   //
-  Int_t    idx = (d == 1 ? 0 : 2*(d - 2) + 1 + ((r=='I' || r=='i') ? 0 : 1));
-  if (fMultCuts[idx] > 0) return fMultCuts[idx];
-  if (fMultCut > 0) return fMultCut;
+  Double_t rcut = GetRingHistos(d, r)->fMultCut;
+  // Int_t    idx = (d == 1 ? 0 : 2*(d - 2) + 1 + ((r=='I' || r=='i') ? 0 : 1));
+  // if (fMultCuts[idx] > 0) return fMultCuts[idx];
+  if (rcut > 0) return rcut;
 
   AliForwardCorrectionManager&  fcm = AliForwardCorrectionManager::Instance();
   AliFMDCorrELossFit* fits = fcm.GetELossFit();
@@ -306,6 +307,10 @@ AliFMDDensityCalculator::GetMultCut(UShort_t d, Char_t r, Double_t eta,
   // Return:
   //    Lower cut on multiplicity
   //
+  Double_t rcut = GetRingHistos(d, r)->fMultCut;
+  if (rcut > 0) return rcut;
+  if (fMultCut > 0) return fMultCut;
+
   AliForwardCorrectionManager&  fcm = AliForwardCorrectionManager::Instance();
   AliFMDCorrELossFit* fits = fcm.GetELossFit();
   Int_t iEta = fits->FindEtaBin(eta);
@@ -367,6 +372,7 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
          if (cut > 0 && mult > cut) 
            n = NParticles(mult,d,r,s,t,vtxbin,eta,lowFlux);
          
+         rh->fELoss->Fill(mult);
          rh->fEvsN->Fill(mult,n);
          rh->fEtaVsN->Fill(eta, n);
          
@@ -376,7 +382,8 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
          rh->fEvsM->Fill(mult,n);
          rh->fEtaVsM->Fill(eta, n);
          rh->fCorr->Fill(eta, c);
-         
+
+         if (n > 0.9 && c > 0) rh->fELossUsed->Fill(mult);
          if (n > 0.9 && c > 0) rh->fBasicHits->Fill(eta,phi, 1./c);
          else                  rh->fEmptyStrips->Fill(eta,phi);
            
@@ -826,6 +833,25 @@ AliFMDDensityCalculator::DefineOutput(TList* dir)
   d->Add(fMaxWeights);
   d->Add(fLowCuts);
 
+  TNamed* sigma  = new TNamed("sigma",
+                              (fIncludeSigma ? "included" : "excluded"));
+  TNamed* maxP   = new TNamed("maxParticle", Form("%d", fMaxParticles));
+  TNamed* method = new TNamed("method", 
+                             (fUsePoisson ? "Poisson" : "Energy loss"));
+  TNamed* phiA   = new TNamed("phiAcceptance", 
+                             (fUsePhiAcceptance ? "enabled" : "disabled"));
+  TNamed* etaL   = new TNamed("etaLumping", Form("%d", fEtaLumping));
+  TNamed* phiL   = new TNamed("phiLumping", Form("%d", fPhiLumping));
+  TParameter<double>* nxi = new TParameter<double>("nXi", fNXi);
+
+  d->Add(sigma);
+  d->Add(maxP);
+  d->Add(method);
+  d->Add(phiA);
+  d->Add(etaL);
+  d->Add(phiL);
+  d->Add(nxi);
+
   TIter    next(&fRingHistos);
   RingHistos* o = 0;
   while ((o = static_cast<RingHistos*>(next()))) {
@@ -904,7 +930,10 @@ AliFMDDensityCalculator::RingHistos::RingHistos()
     fTotalStrips(0),
     fEmptyStrips(0),
     fBasicHits(0),
-    fEmptyVsTotal(0)
+    fEmptyVsTotal(0),
+    fELoss(0),
+    fELossUsed(0),
+    fMultCut(0)
 {
   // 
   // Default CTOR
@@ -924,7 +953,9 @@ AliFMDDensityCalculator::RingHistos::RingHistos(UShort_t d, Char_t r)
     fTotalStrips(0),
     fEmptyStrips(0),
     fBasicHits(0),
-    fEmptyVsTotal(0)
+    fEmptyVsTotal(0),
+    fELossUsed(0),
+    fMultCut(0)
 {
   // 
   // Constructor
@@ -994,6 +1025,24 @@ AliFMDDensityCalculator::RingHistos::RingHistos(UShort_t d, Char_t r)
   fEmptyVsTotal->SetXTitle("Total # strips");
   fEmptyVsTotal->SetYTitle("Empty # strips");
   fEmptyVsTotal->SetZTitle("Correlation");  
+
+  fELoss = new TH1D("eloss", "#Delta/#Delta_{mip} in all strips", 
+                   600, 0, 15);
+  fELoss->SetXTitle("#Delta/#Delta_{mip} (selected)");
+  fELoss->SetYTitle("P(#Delta/#Delta_{mip})");
+  fELoss->SetFillColor(Color()-2);
+  fELoss->SetFillStyle(3003);
+  fELoss->SetLineColor(kBlack);
+  fELoss->SetLineStyle(2);
+  fELoss->SetLineWidth(2);
+  fELoss->SetDirectory(0);
+
+  fELossUsed = static_cast<TH1D*>(fELoss->Clone("elossUsed"));
+  fELossUsed->SetTitle("#Delta/#Delta_{mip} in used strips");
+  fELossUsed->SetFillStyle(3002);
+  fELossUsed->SetLineStyle(1);
+  fELossUsed->SetDirectory(0);
+  
 }
 //____________________________________________________________________
 AliFMDDensityCalculator::RingHistos::RingHistos(const RingHistos& o)
@@ -1008,7 +1057,10 @@ AliFMDDensityCalculator::RingHistos::RingHistos(const RingHistos& o)
     fTotalStrips(o.fTotalStrips),
     fEmptyStrips(o.fEmptyStrips),
     fBasicHits(o.fBasicHits),
-    fEmptyVsTotal(o.fEmptyVsTotal)
+    fEmptyVsTotal(o.fEmptyVsTotal),
+    fELoss(o.fELoss),
+    fELossUsed(o.fELossUsed),
+    fMultCut(o.fMultCut)
 {
   // 
   // Copy constructor 
@@ -1049,9 +1101,11 @@ AliFMDDensityCalculator::RingHistos::operator=(const RingHistos& o)
   fEtaVsM         = static_cast<TProfile*>(o.fEtaVsM->Clone());
   fCorr           = static_cast<TProfile*>(o.fCorr->Clone());
   fDensity        = static_cast<TH2D*>(o.fDensity->Clone());
-  fELossVsPoisson = static_cast<TH2D*>(o.fELossVsPoisson);
-  fTotalStrips    = static_cast<TH2D*>(o.fTotalStrips);
-  fEmptyStrips    = static_cast<TH2D*>(o.fEmptyStrips);
+  fELossVsPoisson = static_cast<TH2D*>(o.fELossVsPoisson->Clone());
+  fTotalStrips    = static_cast<TH2D*>(o.fTotalStrips->Clone());
+  fEmptyStrips    = static_cast<TH2D*>(o.fEmptyStrips->Clone());
+  fELoss          = static_cast<TH1D*>(o.fELoss->Clone());
+  fELossUsed      = static_cast<TH1D*>(o.fELossUsed->Clone());
   
   return *this;
 }
@@ -1061,15 +1115,6 @@ AliFMDDensityCalculator::RingHistos::~RingHistos()
   // 
   // Destructor 
   //
-  if (fEvsN)           delete fEvsN;
-  if (fEvsM)           delete fEvsM;
-  if (fEtaVsN)         delete fEtaVsN;
-  if (fEtaVsM)         delete fEtaVsM;
-  if (fCorr)           delete fCorr;
-  if (fDensity)        delete fDensity;
-  if (fELossVsPoisson) delete fELossVsPoisson;
-  if (fTotalStrips)    delete fTotalStrips;
-  if (fEmptyStrips)    delete fEmptyStrips;
 }
 
 //____________________________________________________________________
@@ -1144,11 +1189,14 @@ AliFMDDensityCalculator::RingHistos::Output(TList* dir)
   d->Add(fDensity);
   d->Add(fELossVsPoisson);
   d->Add(fEmptyVsTotal);
+  d->Add(fELoss);
+  d->Add(fELossUsed);
   // d->Add(fTotalStrips);
   // d->Add(fEmptyStrips);
   // d->Add(fBasicHits);
   
-  
+  TParameter<double>* cut = new TParameter<double>("cut", fMultCut);
+  d->Add(cut);
 }
 
 //____________________________________________________________________
index 80309cd..7624a14 100644 (file)
@@ -368,10 +368,12 @@ protected:
     TH2D*     fTotalStrips;    // Total number of strips in a region
     TH2D*     fEmptyStrips;    // Total number of strips in a region
     TH2D*     fBasicHits  ;    // Total number basic hits in a region
-    TH2D*     fEmptyVsTotal;   // # of empty strips vs total number of strips 
+    TH2D*     fEmptyVsTotal;   // # of empty strips vs total number of # strips 
+    TH1D*     fELoss;          // Energy loss as seen by this 
+    TH1D*     fELossUsed;      // Energy loss in strips with signal 
+    Double_t  fMultCut;        // If set, use this
     
-    
-    ClassDef(RingHistos,4);
+    ClassDef(RingHistos,5);
   };
   /** 
    * Get the ring histogram container 
@@ -405,9 +407,9 @@ protected:
   Int_t    fEtaLumping;    //  How to lump eta bins for Poisson 
   Int_t    fPhiLumping;    //  How to lump phi bins for Poisson 
   Int_t    fDebug;         //  Debug level 
-  Double_t fMultCuts[5];   //  Per-ring multiplicity cuts
+  // Double_t fMultCuts[5];   //  Per-ring multiplicity cuts
 
-  ClassDef(AliFMDDensityCalculator,3); // Calculate Nch density 
+  ClassDef(AliFMDDensityCalculator,4); // Calculate Nch density 
 };
 
 #endif
index adb980c..7606cc1 100644 (file)
@@ -395,6 +395,7 @@ AliFMDEnergyFitter::MakeCorrectionsObject(TList* d)
   TString fileName(mgr.GetFilePath(AliForwardCorrectionManager::kELossFits));
   AliInfo(Form("Object %s created in output - should be extracted and copied "
               "to %s", oName.Data(), fileName.Data()));
+  d->Add(new TNamed("filename", fileName.Data()));
   d->Add(obj, oName.Data());
 }
 
index a3067c7..4a91981 100644 (file)
@@ -145,13 +145,15 @@ AliFMDEnergyFitterTask::InitializeSubs()
   // 
   //
   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
-  fcm.Init(fEventInspector.GetCollisionSystem(), 
-          fEventInspector.GetEnergy(),
-          fEventInspector.GetField(), 0);
+  UShort_t sys = fEventInspector.GetCollisionSystem();
+  UShort_t sNN = fEventInspector.GetEnergy();
+  Short_t  fld = fEventInspector.GetField();
+  fcm.Init(sys, sNN, fld, 0);
   TAxis eAxis(0,0,0);
   TAxis vAxis(10,-10,10);
   fEnergyFitter.Init(eAxis);
   fEventInspector.Init(vAxis);
+
 }
 
 //____________________________________________________________________
index 9a19f36..ced2c90 100644 (file)
@@ -314,6 +314,29 @@ AliFMDEventInspector::Init(const TAxis& vtxAxis)
 
 //____________________________________________________________________
 void
+AliFMDEventInspector::StoreInformation()
+{
+  // Write TNamed objects to output list containing information about
+  // the running conditions 
+  if (!fList) return;
+
+  TNamed* sys = new TNamed("sys", "");
+  TNamed* sNN = new TNamed("sNN", "");
+  TNamed* fld = new TNamed("field", "");
+  sys->SetTitle(AliForwardUtil::CollisionSystemString(fCollisionSystem));
+  sNN->SetTitle(AliForwardUtil::CenterOfMassEnergyString(fEnergy));
+  fld->SetTitle(AliForwardUtil::MagneticFieldString(fField));
+  sys->SetUniqueID(fCollisionSystem);
+  sNN->SetUniqueID(fEnergy);
+  fld->SetUniqueID(fField);
+
+  fList->Add(sys);
+  fList->Add(sNN);
+  fList->Add(fld);
+}
+
+//____________________________________________________________________
+void
 AliFMDEventInspector::DefineOutput(TList* dir)
 {
   // 
@@ -751,7 +774,8 @@ AliFMDEventInspector::ReadRunDetails(const AliESDEvent* esd)
                                            2 * esd->GetBeamEnergy());
   fField           = 
     AliForwardUtil::ParseMagneticField(esd->GetMagneticField());
-  
+
+  StoreInformation();
   if (fCollisionSystem   == AliForwardUtil::kUnknown || 
       fEnergy            <= 0                        || 
       TMath::Abs(fField) >  10) 
index c526202..9730200 100644 (file)
@@ -221,6 +221,19 @@ public:
    * @param option Not used 
    */
   void Print(Option_t* option="") const;
+  /** 
+   * Store information about running conditions in output list 
+   * 
+   * 3 TNamed objects are defined.   The names are fixed, but the
+   * title is a string representation of the information, and the
+   * unique ID contains the identifier 
+   *
+   * - sys   Contains the collision system string and identifier. 
+   * - sNN   Contains the center-of-mass energy per nucleon (GeV)
+   * - field Contains the L3 magnetic field (kG)
+   * 
+   */
+  virtual void StoreInformation();
 protected:
   /** 
    * Read the trigger information from the ESD event 
index 507fcbd..6ced6d8 100644 (file)
@@ -46,7 +46,8 @@ AliFMDMCEventInspector::AliFMDMCEventInspector()
     fHBvsCent(0),
     fHVzComp(0),
     fHCentVsPart(0),
-    fHCentVsBin(0)
+    fHCentVsBin(0),
+    fProduction("")
 {
   // 
   // Constructor 
@@ -64,7 +65,8 @@ AliFMDMCEventInspector::AliFMDMCEventInspector(const char* /* name */)
     fHBvsCent(0),
     fHVzComp(0),
     fHCentVsPart(0),
-    fHCentVsBin(0)
+    fHCentVsBin(0),
+    fProduction("")
 {
   // 
   // Constructor 
@@ -85,7 +87,8 @@ AliFMDMCEventInspector::AliFMDMCEventInspector(const AliFMDMCEventInspector& o)
     fHBvsCent(0),
     fHVzComp(0),
     fHCentVsPart(0),
-    fHCentVsBin(0)
+    fHCentVsBin(0),
+    fProduction("")
 {
   // 
   // Copy constructor 
@@ -218,6 +221,63 @@ AliFMDMCEventInspector::Init(const TAxis& vtxAxis)
 }
 
 //____________________________________________________________________
+void
+AliFMDMCEventInspector::StoreInformation()
+{
+  // Store information about running conditions in the output list 
+  if (!fList) return;
+  AliFMDEventInspector::StoreInformation();
+  TNamed* mc = new TNamed("mc", fProduction.Data());
+  mc->SetUniqueID(1);
+  fList->Add(mc);
+}
+
+namespace
+{
+  TString& AppendField(TString& s, bool test, const char* f)
+  {
+    if (!test) return s;
+    if (!s.IsNull()) s.Append(",");
+    s.Append(f);
+    return s;
+  }
+}
+
+//____________________________________________________________________
+void
+AliFMDMCEventInspector::ReadProductionDetails(AliMCEvent* event)
+{
+  //Assign MC only triggers : True NSD etc.
+  AliHeader*               header          = event->Header();
+  AliGenEventHeader*       genHeader       = header->GenEventHeader();
+  AliCollisionGeometry*    colGeometry     = 
+    dynamic_cast<AliCollisionGeometry*>(genHeader);
+  AliGenPythiaEventHeader* pythiaHeader    = 
+    dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
+  AliGenDPMjetEventHeader* dpmHeader       = 
+    dynamic_cast<AliGenDPMjetEventHeader*>(genHeader);
+  AliGenGeVSimEventHeader* gevHeader       = 
+    dynamic_cast<AliGenGeVSimEventHeader*>(genHeader);
+  AliGenHerwigEventHeader* herwigHeader    = 
+    dynamic_cast<AliGenHerwigEventHeader*>(genHeader);
+  AliGenHijingEventHeader* hijingHeader    = 
+    dynamic_cast<AliGenHijingEventHeader*>(genHeader);
+  // AliGenHydjetEventHeader* hydjetHeader    = 
+  //   dynamic_cast<AliGenHydjetEventHeader*>(genHeader);
+  AliGenEposEventHeader*   eposHeader      = 
+    dynamic_cast<AliGenEposEventHeader*>(genHeader);
+
+  AppendField(fProduction, colGeometry,  "Geometry");
+  AppendField(fProduction, pythiaHeader, "Pythia");
+  AppendField(fProduction, dpmHeader,    "DPM");
+  AppendField(fProduction, gevHeader,    "GevSim");
+  AppendField(fProduction, herwigHeader, "Herwig");
+  AppendField(fProduction, hijingHeader, "Hijing");
+  // AppendField(fProduction, hydjetHeader, "Hydjet");
+  AppendField(fProduction, eposHeader,   "EPOS");
+}
+
+//____________________________________________________________________
 UInt_t
 AliFMDMCEventInspector::ProcessMC(AliMCEvent*       event, 
                                  UInt_t&           triggers,
index f84cd96..c339c8c 100644 (file)
@@ -117,6 +117,22 @@ public:
   virtual Bool_t CompareResults(Double_t vz,    Double_t trueVz, 
                                Double_t cent,  Double_t b,
                                Int_t    npart, Int_t    nbin);
+  /** 
+   * Store information about running conditions in output list 
+   * 
+   * The 3 TNamed objects from AliFMDEventInspector::StoreInformation
+   * are defined.  In addition, a fourth TNamed object is defined.
+   * The presence of this indicate MC data.
+   *
+   * - mc    Nothing special, and unique id is 1
+   */
+  virtual void StoreInformation();
+  /** 
+   * Read the production details 
+   * 
+   * @param event MC event
+   */
+  virtual void ReadProductionDetails(AliMCEvent* event);
 protected:
   /** 
    * Read centrality from event 
@@ -150,7 +166,8 @@ protected:
   TH2F* fHVzComp;  // True vs reconstructed vz
   TH2F* fHCentVsPart; // Centrality versus # participants 
   TH2F* fHCentVsBin;  // Centrality versus # binary collisions 
-  ClassDef(AliFMDMCEventInspector,2); // Inspect the event 
+  TString fProduction; // Production information 
+  ClassDef(AliFMDMCEventInspector,3); // Inspect the event 
 };
 
 #endif
index c6feb6d..66e5362 100644 (file)
@@ -32,6 +32,7 @@
 #include <AliLog.h>
 #include <TROOT.h>
 #include <THStack.h>
+#include <TParameter.h>
 #include <iostream>
 #include <iomanip>
 
@@ -59,7 +60,7 @@ AliFMDSharingFilter::AliFMDSharingFilter()
     fHighCuts(0),
     fOper(0),
     fDebug(0),
-    fZeroSharedHitsBelowThreshold(0)
+    fZeroSharedHitsBelowThreshold(false)
 {
   // 
   // Default Constructor - do not use 
@@ -78,7 +79,7 @@ AliFMDSharingFilter::AliFMDSharingFilter(const char* title)
     fHighCuts(0),
     fOper(0),
     fDebug(0),
-    fZeroSharedHitsBelowThreshold(0)
+    fZeroSharedHitsBelowThreshold(false)
 {
   // 
   // Constructor 
@@ -558,9 +559,7 @@ AliFMDSharingFilter::MultiplicityOfStrip(Double_t thisE,
   DBGL(3, Form(" %9f<%9f<%9f (was %9f), zeroed, candidate", 
                   lowCut, totalE, highCut, thisE));
   thisStatus = kCandidate;
-  if(fZeroSharedHitsBelowThreshold)
-    return 0;
-  else return totalE; 
+  return (fZeroSharedHitsBelowThreshold ? 0 : totalE); 
 }
           
 //_____________________________________________________________________
@@ -784,6 +783,25 @@ AliFMDSharingFilter::DefineOutput(TList* dir)
   fHighCuts->SetDirectory(0);
   d->Add(fHighCuts);
 
+  TParameter<double>* lowCut = new TParameter<double>("lowCut", fLowCut);
+  TParameter<double>* nXi    = new TParameter<double>("nXi", fNXi);
+  TNamed*             sigma  = new TNamed("sigma", fIncludeSigma ? 
+                                         "included" : "excluded");
+  sigma->SetUniqueID(fIncludeSigma);
+  TNamed*             angle  = new TNamed("angle", fCorrectAngles ? 
+                                         "corrected" : "uncorrected");
+  angle->SetUniqueID(fCorrectAngles);
+  TNamed*             low    = new TNamed("lowSignal", 
+                                         fZeroSharedHitsBelowThreshold ? 
+                                         "zeroed" : "kept");
+  low->SetUniqueID(fZeroSharedHitsBelowThreshold);
+  d->Add(lowCut);
+  d->Add(nXi);
+  d->Add(sigma);
+  d->Add(angle);
+  d->Add(low);
+  
+
 
   TIter    next(&fRingHistos);
   RingHistos* o = 0;
@@ -810,7 +828,9 @@ AliFMDSharingFilter::Print(Option_t* /*option*/) const
            << ind << " Low cut:                " << fLowCut << '\n'
            << ind << " N xi factor:            " << fNXi    << '\n'
            << ind << " Include sigma in cut:   " << fIncludeSigma << '\n'
-           << ind << " Use corrected angles:   " << fCorrectAngles 
+           << ind << " Use corrected angles:   " << fCorrectAngles << '\n'
+           << ind << " Zero below threshold:   " 
+           << fZeroSharedHitsBelowThreshold 
            << std::noboolalpha << std::endl;
 }
   
@@ -852,16 +872,19 @@ AliFMDSharingFilter::RingHistos::RingHistos(UShort_t d, Char_t r)
   //    r ring 
   //
   fBefore = new TH1D("esdEloss", "Energy loss (reconstruction)", 
-                    1000, 0, 25);
+                    600, 0, 15);
   fBefore->SetXTitle("#Delta E/#Delta E_{mip}");
   fBefore->SetYTitle("P(#Delta E/#Delta E_{mip})");
   fBefore->SetFillColor(Color());
   fBefore->SetFillStyle(3001);
+  fBefore->SetLineColor(kBlack);
+  fBefore->SetLineStyle(2);
   fBefore->SetDirectory(0);
 
   fAfter  = static_cast<TH1D*>(fBefore->Clone("anaEloss"));
   fAfter->SetTitle("Energy loss in %s (sharing corrected)");
   fAfter->SetFillColor(Color()+2);
+  fAfter->SetLineStyle(1);
   fAfter->SetDirectory(0);
 
   Double_t max = 15;
@@ -875,7 +898,7 @@ AliFMDSharingFilter::RingHistos::RingHistos(UShort_t d, Char_t r)
   fBeforeAfter->SetDirectory(0);
 
   fNeighborsBefore = static_cast<TH2D*>(fBeforeAfter->Clone("neighborsBefore"));
-  fNeighborsBefore->SetTitle("Correlation of neighbors befire");
+  fNeighborsBefore->SetTitle("Correlation of neighbors before");
   fNeighborsBefore->SetXTitle("#Delta E_{i}/#Delta E_{mip}");
   fNeighborsBefore->SetYTitle("#Delta E_{i+1}/#Delta E_{mip}");
   fNeighborsBefore->SetDirectory(0);
@@ -984,11 +1007,14 @@ AliFMDSharingFilter::RingHistos::ScaleHistograms(const TList* dir,
   TList* l = GetOutputList(dir);
   if (!l) return; 
 
-  TH1D* before = static_cast<TH1D*>(l->FindObject("esdELoss"));
-  TH1D* after  = static_cast<TH1D*>(l->FindObject("anaELoss"));
-  TH2D* summed = static_cast<TH2D*>(l->FindObject("summed"));
+#if 0
+  TH1D* before = static_cast<TH1D*>(l->FindObject("esdEloss"));
+  TH1D* after  = static_cast<TH1D*>(l->FindObject("anaEloss"));
   if (before) before->Scale(1./nEvents);
   if (after)  after->Scale(1./nEvents);
+#endif
+
+  TH2D* summed = static_cast<TH2D*>(l->FindObject("summed"));
   if (summed) summed->Scale(1./nEvents);
   
 }
index cdcaadf..0ea7f56 100644 (file)
@@ -289,6 +289,10 @@ AliForwardMCMultiplicityTask::UserExec(Option_t*)
   //    option Not used
   //  
 
+  // Read production details 
+  if (fFirstEvent) 
+    fEventInspector.ReadProductionDetails(MCEvent());
+    
   // Get the input data 
   AliESDEvent* esd     = GetESDEvent();
   AliMCEvent*  mcEvent = MCEvent();
index 9b44572..7cd4e80 100644 (file)
@@ -279,8 +279,8 @@ AliForwardMultiplicityTask::UserExec(Option_t*)
   }
 
   // Calculate the inclusive charged particle density 
-  //if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) { 
-  if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) { 
+  if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) { 
+    // if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) { 
     AliWarning("Density calculator failed!");
     return;
   }
index 490f7b2..93e05ca 100644 (file)
@@ -217,13 +217,14 @@ AliMCTruthdNdetaTask::CentralityBin::End(TList*      sums,
                                        shapeCorr, trigEff, 
                                        symmetrice, rebin, 
                                        rootProj, corrEmpty, cutEdges,
-                                       triggerMask, marker, color, mclist, truthlist);
+                                       triggerMask, marker, color, mclist, 
+                                       truthlist);
 
   fSumTruth     = static_cast<TH2D*>(fSums->FindObject("truth"));
   
 
   if (fSumTruth) { 
-    Int_t n0 = fSumTruth->GetBinContent(0,0);
+    Int_t n0 = Int_t(fSumTruth->GetBinContent(0,0));
     Int_t n  = (triggerMask == AliAODForwardMult::kNSD ? 
                Int_t(fTriggers->GetBinContent(AliAODForwardMult::kBinMCNSD)) : 
                Int_t(fTriggers->GetBinContent(AliAODForwardMult::kBinAll)));
index b917273..31d6600 100644 (file)
@@ -43,8 +43,8 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
     mcTask->SetOnlyPrimary(true);
   }
 #endif
-  Double_t nXi = 1; // mc ? 1 : .5;
-  Bool_t   includeSigma = true;
+  Double_t nXi = 2; // mc ? 1 : .5;
+  Bool_t   includeSigma = false; 
 
   // --- Event inspector ---------------------------------------------
   // Set the number of SPD tracklets for which we consider the event a
@@ -55,13 +55,14 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
 
   // --- Sharing filter ----------------------------------------------
   // Set the low cut used for sharing - overrides settings in eloss fits
-  task->GetSharingFilter().SetLowCut(mc ? 0.15 : 0.2);
+  task->GetSharingFilter().SetLowCut(0.15); // mc ? 0.15 : 0.2);
   // Set the number of xi's (width of landau peak) to stop at 
   task->GetSharingFilter().SetNXi(nXi);
   // Set whether or not to include sigma in cut
   task->GetSharingFilter().SetIncludeSigma(includeSigma);
   // Enable use of angle corrected signals in the algorithm 
   task->GetSharingFilter().SetUseAngleCorrectedSignals(true);
+  task->GetSharingFilter().SetZeroSharedHitsBelowThreshold(false);
 
   // --- Density calculator 
   // Set the maximum number of particle to try to reconstruct 
@@ -69,9 +70,10 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
   // Wet whether to use poisson statistics to estimate N_ch
   task->GetDensityCalculator().SetUsePoisson(true);
   // Set the lower multiplicity cut.  Overrides setting in energy loss fits.
-  task->GetDensityCalculator().SetMultCut(.3); //was 0.3
+  task->GetDensityCalculator().SetMultCut(-1); //was 0.3
   // Set the lower per-ring multiplicity cuts 
-  task->GetDensityCalculator().SetMultCuts(-1,-1,-1,-1,-1);
+  task->GetDensityCalculator().SetMultCuts(0.3,0.3,0.3,0.3,0.3);
+  // task->GetDensityCalculator().SetMultCuts(-1,-1,-1,-1,-1);
   // USe this many times xi+sigma below MPV 
   task->GetDensityCalculator().SetNXi(nXi);
   // Set whether or not to include sigma in cut
index a23ebfd..0974543 100644 (file)
@@ -116,7 +116,7 @@ void MakeAOD(const char* esddir,
     static_cast<AliPhysicsSelection*>(ih->GetEventSelection());
   // Ignore trigger class when selecting events.  This mean that we
   // get offline+(A,C,E) events too
-  ps->SetSkipTriggerClassSelection(true);
+  // ps->SetSkipTriggerClassSelection(true);
   
 
 #if 0
index 2bbbd80..2abf54d 100644 (file)
@@ -1238,10 +1238,29 @@ GetSingle(UShort_t which,
       }
     }
   }
+#if 0
+  if (ret) {
+    if (!wn.IsNull()) wn.Append(",");
+    switch (which) { 
+    case PYTHIA: wn.Append("Pythia"); break;
+    case UA5:    wn.Append("UA5");    break;
+    case CMS:    wn.Append("CMS");    break;
+    case ALICE:  wn.Append("ALICE");  break;
+    }      
+  }
+#endif
   return ret;
 }
 
-         
+//____________________________________________________________________
+TString&
+AppendItem(TString& s, char delim, const char* what)     
+{
+  if (!s.IsNull()) s.Append(Form("%c", delim));
+  s.Append(what);
+  return s;
+}
+
 //____________________________________________________________________
 /** 
  * Get a multi graph of data for a given energy and trigger type 
@@ -1293,6 +1312,10 @@ GetData(UShort_t sys,
     if (!(type & 0x7)) 
       Warning("GetData", "Unknown trigger mask 0x%x", type);
 
+    if (TMath::Abs(energy-2750) < 11) {
+      Warning("GetData", "Using 2360GeV data for %dGeV comparison", energy);
+      energy = 2360;
+    }
     if (!(TMath::Abs(energy-900) < 10 || 
          TMath::Abs(energy-2360) < 10 || 
          TMath::Abs(energy-7000) < 10)) {
@@ -1303,19 +1326,19 @@ GetData(UShort_t sys,
     
     sn = "pp";
 
-    if (type & 0x1) tn.Append("INEL");
-    if (type & 0x2) { if (!tn.IsNull()) tn.Append("|"); tn.Append("INEL>0"); }
-    if (type & 0x4) { if (!tn.IsNull()) tn.Append("|"); tn.Append("NSD"); }
+    if (type & 0x1) AppendItem(tn, '|', "INEL");
+    if (type & 0x2) AppendItem(tn, '|', "INEL>0");
+    if (type & 0x4) AppendItem(tn, '|', "NSD");
 
     Bool_t seenUA5 = false;
     for (Int_t i = 0; i < 3; i++) { 
       UShort_t mask = (1 << i);
       if ((type & mask) == 0) continue;
-      TGraphAsymmErrors* gUAp = (ua5   ? GetSingle(UA5,   sys,energy,mask): 0);
-      TGraphAsymmErrors* gUAn = (ua5   ? GetSingle(UA5+10,sys,energy,mask): 0);
-      TGraphAsymmErrors* gCMS = (cms   ? GetSingle(CMS,   sys,energy,mask): 0);
-      TGraphAsymmErrors* gALI = (alice ? GetSingle(ALICE, sys,energy,mask): 0);
-      TGraphAsymmErrors* gPYT = (pythia? GetSingle(PYTHIA,sys,energy,mask): 0);
+      TGraphAsymmErrors* gUAp =(ua5   ?GetSingle(UA5,   sys,energy,mask):0);
+      TGraphAsymmErrors* gUAn =(ua5   ?GetSingle(UA5+10,sys,energy,mask):0);
+      TGraphAsymmErrors* gCMS =(cms   ?GetSingle(CMS,   sys,energy,mask):0);
+      TGraphAsymmErrors* gALI =(alice ?GetSingle(ALICE, sys,energy,mask):0);
+      TGraphAsymmErrors* gPYT =(pythia?GetSingle(PYTHIA,sys,energy,mask):0);
       if (gUAp) mp->Add(gUAp);
       if (gUAn) mp->Add(gUAn);
       if (gCMS) mp->Add(gCMS);
index da4ab72..00fc109 100644 (file)
@@ -1761,7 +1761,7 @@ protected:
     LoadLibrary("PWG2forward2", mode, par, true);
     
     // --- Set load path ---------------------------------------------
-    gROOT->SetMacroPath(Form(".:%s:$(ALICE_ROOT)/PWG2/FORWARD/analysis2",
+    gROOT->SetMacroPath(Form("%s:$(ALICE_ROOT)/PWG2/FORWARD/analysis2",
                             gROOT->GetMacroPath()));
 
     // --- Check if this is MC ---------------------------------------
@@ -1799,11 +1799,11 @@ protected:
     AliPhysicsSelection* ps = 
       dynamic_cast<AliPhysicsSelection*>(ih->GetEventSelection());
     if (!ps) 
-      Fatal("CreatePhysicsSelection", "Couldn't get PhysicsSelection (%p)", ps);
+      Fatal("CreatePhysicsSelection", "Couldn't get PhysicsSelection (%p)",ps);
 
     // --- Ignore trigger class when selecting events.  This means ---
     // --- that we get offline+(A,C,E) events too --------------------
-    ps->SetSkipTriggerClassSelection(true);
+    // ps->SetSkipTriggerClassSelection(true);
   }
   //__________________________________________________________________
   /** 
@@ -2071,7 +2071,7 @@ protected:
     LoadLibrary("PWG2forward2", mode, par, true);
     
     // --- Set load path ---------------------------------------------
-    gROOT->SetMacroPath(Form(".:%s:$(ALICE_ROOT)/PWG2/FORWARD/analysis2",
+    gROOT->SetMacroPath(Form("%s:$(ALICE_ROOT)/PWG2/FORWARD/analysis2",
                             gROOT->GetMacroPath()));
 
     // --- Check if this is MC ---------------------------------------