Added a prototype of V0 cut.
authorpulvir <pulvir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Jun 2011 16:16:25 +0000 (16:16 +0000)
committerpulvir <pulvir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Jun 2011 16:16:25 +0000 (16:16 +0000)
Updated the example macros

13 files changed:
PWG2/CMakelibPWG2resonances.pkg
PWG2/PWG2resonancesLinkDef.h
PWG2/RESONANCES/AliRsnCutV0.cxx [new file with mode: 0644]
PWG2/RESONANCES/AliRsnCutV0.h [new file with mode: 0644]
PWG2/RESONANCES/AliRsnTarget.h
PWG2/RESONANCES/macros/mini/AddAnalysisTaskRsnMiniTest.C [new file with mode: 0644]
PWG2/RESONANCES/macros/mini/AnalysisSetupRsnMini.C
PWG2/RESONANCES/macros/mini/ConfigKStarSimple.C [new file with mode: 0644]
PWG2/RESONANCES/macros/mini/ConfigPhiSimple.C [new file with mode: 0644]
PWG2/RESONANCES/macros/mini/ConfigSigmaStar.C [new file with mode: 0644]
PWG2/RESONANCES/macros/mini/ConfigSigmaStarMC.C [new file with mode: 0644]
PWG2/RESONANCES/macros/mini/runLocal.C
PWG2/RESONANCES/macros/mini/runPluginProof.C

index 9ecb595..dc022c8 100644 (file)
@@ -51,6 +51,7 @@ set ( SRCS RESONANCES/AliRsnDaughter.cxx
            RESONANCES/AliRsnCutPion2010PP.cxx
            RESONANCES/AliRsnCutProton2010PP.cxx
            RESONANCES/AliRsnCutDaughterKStar2010PP.cxx 
+           RESONANCES/AliRsnCutV0.cxx 
            RESONANCES/AliRsnCutSet.cxx
            RESONANCES/AliRsnExpression.cxx
            RESONANCES/AliRsnVariableExpression.cxx
index b72368f..15e45e1 100644 (file)
@@ -30,7 +30,7 @@
 #pragma link C++ class AliRsnCutPion2010PP+;
 #pragma link C++ class AliRsnCutProton2010PP+;
 #pragma link C++ class AliRsnCutDaughterKStar2010PP+;
-
+#pragma link C++ class AliRsnCutV0+;
 
 #pragma link C++ class AliRsnCutSet+;
 #pragma link C++ class AliRsnExpression+;
diff --git a/PWG2/RESONANCES/AliRsnCutV0.cxx b/PWG2/RESONANCES/AliRsnCutV0.cxx
new file mode 100644 (file)
index 0000000..3cc1e8b
--- /dev/null
@@ -0,0 +1,209 @@
+//
+// Class AliRsnCutV0
+//
+// General implementation of a single cut strategy, which can be:
+// - a value contained in a given interval  [--> IsBetween()   ]
+// - a value equal to a given reference     [--> MatchesValue()]
+//
+// In all cases, the reference value(s) is (are) given as data members
+// and each kind of cut requires a given value type (Int, UInt, Double),
+// but the cut check procedure is then automatized and chosen thanks to
+// an enumeration of the implemented cut types.
+// At the end, the user (or any other point which uses this object) has
+// to use the method IsSelected() to check if this cut has been passed.
+//
+// authors: Martin Vala (martin.vala@cern.ch)
+//          Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
+//
+
+#include <Riostream.h>
+#include <TFormula.h>
+#include <TBits.h>
+
+#include "AliLog.h"
+#include "AliESDtrackCuts.h"
+
+#include "AliRsnEvent.h"
+#include "AliRsnDaughter.h"
+#include "AliRsnCutV0.h"
+
+ClassImp(AliRsnCutV0)
+
+//_________________________________________________________________________________________________
+AliRsnCutV0::AliRsnCutV0(const char *name, Int_t hypothesis) :
+   AliRsnCut(name, AliRsnTarget::kDaughter),
+   fHypothesis(0),
+   fMass(0.0),
+   fTolerance(0.2),
+   fMaxDCAVertex(0.3),
+   fMinCosPointAngle(0.95),
+   fMaxDaughtersDCA(0.1),
+   fESDtrackCuts(0x0)
+{
+//
+// Default constructor.
+// Initializes all cuts in such a way that all of them are disabled.
+//
+
+   SetHypothesis(hypothesis);
+}
+
+//_________________________________________________________________________________________________
+AliRsnCutV0::AliRsnCutV0(const AliRsnCutV0 &copy) :
+   AliRsnCut(copy),
+   fHypothesis(copy.fHypothesis),
+   fMass(copy.fMass),
+   fTolerance(copy.fTolerance),
+   fMaxDCAVertex(copy.fMaxDCAVertex),
+   fMinCosPointAngle(copy.fMinCosPointAngle),
+   fMaxDaughtersDCA(copy.fMaxDaughtersDCA),
+   fESDtrackCuts(copy.fESDtrackCuts)
+{
+//
+// Copy constructor.
+// Just copy all data member values.
+//
+}
+
+//_________________________________________________________________________________________________
+AliRsnCutV0& AliRsnCutV0::operator=(const AliRsnCutV0 &copy)
+{
+//
+// Assignment operator.
+// Just copy all data member values.
+//
+
+
+   fHypothesis = copy.fHypothesis;
+   fMass = copy.fMass;
+   fTolerance = copy.fTolerance;
+   fMaxDCAVertex = copy.fMaxDCAVertex;
+   fMinCosPointAngle = copy.fMinCosPointAngle;
+   fMaxDaughtersDCA = copy.fMaxDaughtersDCA;
+   fESDtrackCuts = copy.fESDtrackCuts;
+
+   return (*this);
+}
+
+//_________________________________________________________________________________________________
+Bool_t AliRsnCutV0::IsSelected(TObject *object)
+{
+//
+// Cut checker.
+// Checks the type of object being evaluated
+// and then calls the appropriate sub-function (for ESD or AOD)
+//
+
+   // coherence check
+   if (!TargetOK(object)) return kFALSE;
+   
+   // check cast
+   AliESDv0 *v0esd = fDaughter->Ref2ESDv0();
+   AliAODv0 *v0aod = fDaughter->Ref2AODv0();
+   //cout << fDaughter->GetRef()->ClassName() << ' ' << v0esd << ' ' << v0aod << endl;
+   
+   // operate depending on cast
+   if (v0esd) {
+      return CheckESD(v0esd);
+   } else if (v0aod) {
+      return CheckAOD(v0aod);
+   } else {
+      AliDebugClass(1, "Object is not a V0");
+      return kFALSE;
+   }
+}
+
+//_________________________________________________________________________________________________
+Bool_t AliRsnCutV0::CheckESD(AliESDv0 *v0)
+{
+//
+// Check an ESD V0.
+// This is done using the default track checker for ESD.
+// It is declared static, not to recreate it every time.
+//
+
+   AliDebugClass(1, "Check ESD");
+   if (v0->GetOnFlyStatus()) {
+      AliDebugClass(1, "Rejecting V0 in 'on fly' status");
+      return kFALSE; // if kTRUE, then this V0 is recontructed
+   }
+   
+   // retrieve pointer to owner event
+   AliESDEvent *lESDEvent = fEvent->GetRefESD();
+   Double_t xPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetX();
+   Double_t yPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetY();
+   Double_t zPrimaryVertex = lESDEvent->GetPrimaryVertex()->GetZ();
+   AliDebugClass(2, Form("Primary vertex: %f %f %f", xPrimaryVertex, yPrimaryVertex, zPrimaryVertex));
+   
+   // retrieve the V0 daughters
+   UInt_t lIdxPos      = (UInt_t) TMath::Abs(v0->GetPindex());
+   UInt_t lIdxNeg      = (UInt_t) TMath::Abs(v0->GetNindex());
+   AliESDtrack *pTrack = lESDEvent->GetTrack(lIdxPos);
+   AliESDtrack *nTrack = lESDEvent->GetTrack(lIdxNeg);
+   
+   // check quality cuts
+   if (fESDtrackCuts) {
+      AliDebugClass(2, "Checking quality cuts");
+      if (!fESDtrackCuts->IsSelected(pTrack)) {
+         AliDebugClass(2, "Positive daughter failed quality cuts");
+         return kFALSE;
+      }
+      if (!fESDtrackCuts->IsSelected(nTrack)) {
+         AliDebugClass(2, "Negative daughter failed quality cuts");
+         return kFALSE;
+      }
+   }
+   
+   // filter like-sign V0
+   //if ((TMath::Abs(pTrack->GetSign()) - TMath::Abs(nTrack->GetSign()) ) < 0.1) {
+   //   AliDebugClass(2, "Failed like-sign V0 check");
+   //   return kFALSE;
+   //}
+   
+   // check compatibility with expected species hypothesis
+   v0->ChangeMassHypothesis(fHypothesis);
+   if ((TMath::Abs(v0->GetEffMass() - fMass)) > fTolerance) {
+      AliDebugClass(2, "V0 is not in the expected inv mass range");
+      return kFALSE;
+   }
+   
+   // topological checks
+   if (TMath::Abs(v0->GetD(xPrimaryVertex, yPrimaryVertex, zPrimaryVertex)) > fMaxDCAVertex) {
+      AliDebugClass(2, "Failed check on DCA to primary vertes");
+      return kFALSE;
+   }
+   if (TMath::Abs(v0->GetV0CosineOfPointingAngle()) < fMinCosPointAngle) {
+      AliDebugClass(2, "Failed check on cosine of pointing angle");
+      return kFALSE;
+   }
+   if (TMath::Abs(v0->GetDcaV0Daughters()) > fMaxDaughtersDCA) {
+      AliDebugClass(2, "Failed check on DCA between daughters");
+      return kFALSE;
+   }
+   
+   // if we reach this point, all checks were successful
+   AliDebugClass(2, "Good V0 (hallelujah)");
+   return kTRUE;   
+}
+
+//_________________________________________________________________________________________________
+Bool_t AliRsnCutV0::CheckAOD(AliAODv0 *v0)
+{
+//
+// Check an AOD V0.
+// This is done doing directly all checks, since there is not
+// an equivalend checker for AOD tracks
+//
+
+   AliWarning("Cuts is not yet implemented for AOD");
+   
+   return kTRUE;
+}
+
+//_________________________________________________________________________________________________
+void AliRsnCutV0::Print(const Option_t *) const
+{
+//
+// Print information on this cut
+//
+}
diff --git a/PWG2/RESONANCES/AliRsnCutV0.h b/PWG2/RESONANCES/AliRsnCutV0.h
new file mode 100644 (file)
index 0000000..8441cbf
--- /dev/null
@@ -0,0 +1,74 @@
+#ifndef ALIRSNCUTV0_H
+#define ALIRSNCUTV0_H
+
+#include <TMath.h>
+#include <TString.h>
+
+#include "AliRsnCut.h"
+
+class AliESDtrack;
+class AliAODTrack;
+
+class AliRsnCutV0 : public AliRsnCut {
+
+public:
+
+   AliRsnCutV0(const char *name = "AliRsnCutV0", Int_t hypothesis = kLambda0);
+   AliRsnCutV0(const AliRsnCutV0& copy);
+   AliRsnCutV0& operator=(const AliRsnCutV0& copy);
+   virtual ~AliRsnCutV0() { }
+   
+   void           SetESDtrackCuts(AliESDtrackCuts *cuts)   {fESDtrackCuts = cuts;}
+   void           SetHypothesis(Int_t code);
+   void           SetTolerance(Double_t value)             {fTolerance = value;}
+   void           SetMaxDCAVertex(Double_t value)          {fMaxDCAVertex = value;}
+   void           SetMinCosPointingAngle(Double_t value)   {fMinCosPointAngle = value;}
+   void           SetMaxDaughtersDCA(Double_t value)       {fMaxDaughtersDCA = value;}
+
+   virtual Bool_t IsSelected(TObject *obj);
+   virtual void   Print(const Option_t *option = "") const;
+
+protected:
+
+   Bool_t      CheckESD(AliESDv0 *track);
+   Bool_t      CheckAOD(AliAODv0 *track);
+   
+   Int_t            fHypothesis;       // PDG code corresponding to expected V0 hypothesis
+   Double_t         fMass;             // mass corresponding to hypothesis
+   Double_t         fTolerance;        // tolerance in the difference between computed and expected mass
+   Double_t         fMaxDCAVertex;     // max allowed DCA from primary vertex
+   Double_t         fMinCosPointAngle; // min allowed cosine of pointing angle
+   Double_t         fMaxDaughtersDCA;  // max allowed DCA between the two daughers
+   //AliPID::EParticleType fRefPID[2];
+   //Double_t              fNSigmaTPC[2];
+   //Double_t              fNSigmaTOF[3];
+
+   AliESDtrackCuts *fESDtrackCuts;     // quality cuts for v0 daughters
+
+   ClassDef(AliRsnCutV0, 1)
+};
+
+//__________________________________________________________________________________________________
+inline void AliRsnCutV0::SetHypothesis(Int_t code)
+{
+//
+// Assign a V0 species hypothesis, which also assign the expected mass
+//
+
+   fHypothesis = code;
+   
+   switch (fHypothesis) {
+      case kLambda0:
+      case kLambda0Bar:
+         fMass = 1.11568;
+         break;
+      case kK0Short:
+         fMass = 0.497614;
+         break;
+      default:
+         AliError(Form("You are setting an unexpected hypothesis: %d", fHypothesis));
+         fMass = 1E20;
+   }
+}
+
+#endif
index fd188ec..7ec3321 100644 (file)
@@ -12,9 +12,9 @@
 
 #include "TNamed.h"
 
+#include "AliRsnDaughter.h"
 #include "AliRsnEvent.h"
 
-class AliRsnDaughter;
 class AliRsnMother;
 
 class AliRsnTarget : public TNamed {
diff --git a/PWG2/RESONANCES/macros/mini/AddAnalysisTaskRsnMiniTest.C b/PWG2/RESONANCES/macros/mini/AddAnalysisTaskRsnMiniTest.C
new file mode 100644 (file)
index 0000000..1b6e507
--- /dev/null
@@ -0,0 +1,133 @@
+//
+// General macro to configure the RSN analysis task.
+// It calls all configs desired by the user, by means
+// of the boolean switches defined in the first lines.
+// ---
+// Inputs:
+//  1) flag to know if running on MC or data
+//  2) path where all configs are stored
+// ---
+// Returns:
+//  kTRUE  --> initialization successful
+//  kFALSE --> initialization failed (some config gave errors)
+//
+
+Bool_t usePhi   = 1;
+Bool_t useKStar = 1;
+
+AliRsnMiniAnalysisTask * AddAnalysisTaskRsnMiniTest
+(
+   Bool_t      isMC,
+   Bool_t      isPP,
+   const char *path,
+   Int_t       nmix = 0
+)
+{  
+   //
+   // -- INITIALIZATION ----------------------------------------------------------------------------
+   //
+   
+   // retrieve analysis manager
+   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+
+   // create the task and connect with physics selection
+   AliRsnMiniAnalysisTask *task = new AliRsnMiniAnalysisTask("RSN", isMC);
+   mgr->AddTask(task);
+   
+   // settings
+   if (isPP) 
+      task->UseMultiplicity("QUALITY");
+   else
+      task->UseCentrality("V0M");
+   
+   // set mixing
+   task->UseContinuousMix();
+   //task->UseBinnedMix();
+   task->SetNMix(nmix);
+   task->SetMaxDiffVz(1.0);
+   task->SetMaxDiffMult(10.0);
+   task->SetMaxDiffAngle(1E20);
+   
+   //
+   // -- EVENT CUTS (same for all configs) ---------------------------------------------------------
+   //
+   
+   // cut on primary vertex:
+   // - 2nd argument --> |Vz| range
+   // - 3rd argument --> minimum required number of contributors
+   // - 4th argument --> tells if TPC stand-alone vertexes must be accepted
+   AliRsnCutPrimaryVertex *cutVertex = new AliRsnCutPrimaryVertex("cutVertex", 10.0, 0, kFALSE);
+   
+   // set the check for pileup
+   if (isPP) cutVertex->SetCheckPileUp(kTRUE);
+      
+   // define and fill cut set
+   AliRsnCutSet *eventCuts = new AliRsnCutSet("eventCuts", AliRsnTarget::kEvent);
+   eventCuts->AddCut(cutVertex);
+   eventCuts->SetCutScheme(cutVertex->GetName());
+   
+   // set cuts in task
+   task->SetEventCuts(eventCuts);
+   
+   //
+   // -- EVENT-ONLY COMPUTATIONS -------------------------------------------------------------------
+   //
+   
+   // initialize value computation for multiplicity/centrality
+   // second argument tells if the value must be taken from MC
+   // (when this can be done)
+   // after creating the value, the task returns its ID
+   Int_t multID = task->CreateValue(AliRsnMiniValue::kMult, kFALSE);
+   
+   // create event-related output
+   AliRsnMiniOutput *outMult = task->CreateOutput("eventMult", "HIST", "EVENT");
+   // set axes, by passing value ID and defining the binning
+   if (isPP) 
+      outMult->AddAxis(multID, 300, 0.0, 300.0);
+   else
+      outMult->AddAxis(multID, 100, 0.0, 100.0);
+   
+   //
+   // -- PAIR CUTS (common to all resonances) ------------------------------------------------------
+   //
+   
+   AliRsnCutMiniPair *cutY = new AliRsnCutMiniPair("cutRapidity", AliRsnCutMiniPair::kRapidityRange);
+   cutY->SetRangeD(-0.5, 0.5);
+   
+   AliRsnCutSet *cutsPair = new AliRsnCutSet("pairCuts", AliRsnTarget::kMother);
+   cutsPair->AddCut(cutY);
+   cutsPair->SetCutScheme(cutY->GetName());
+   
+   //
+   // -- CONFIGS -----------------------------------------------------------------------------------
+   //
+   
+   if (usePhi) {
+      gROOT->LoadMacro(Form("%s/ConfigPhiSimple.C", path));
+      ConfigPhiSimple(task, isMC, "phi_unlike", "HIST", "PAIR", '+', '-', kTRUE, cutsPair);
+      if (isMC) {
+         gROOT->LoadMacro(Form("%s/ConfigPhiMC.C", path));
+         ConfigPhiMC(task, isPP, "", cutsPair);
+      }
+   }
+   
+   if (useKStar) {
+      gROOT->LoadMacro(Form("%s/ConfigKStarSimple.C", path));
+      ConfigKStarSimple(task, isMC, "kstar_unlike", "HIST", "PAIR", '+', '-', kTRUE, cutsPair);
+      if (isMC) {
+         gROOT->LoadMacro(Form("%s/ConfigKStarMC.C", path));
+         ConfigKStarMC(task, isPP, "", cutsPair);
+      }
+   }
+   
+   //
+   // -- CONTAINERS --------------------------------------------------------------------------------
+   //
+   
+   const char *file = AliAnalysisManager::GetCommonFileName();
+   AliAnalysisDataContainer *output = mgr->CreateContainer("RsnOut", TList::Class(), AliAnalysisManager::kOutputContainer, file);
+   mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
+   mgr->ConnectOutput(task, 1, output);
+
+   return task;
+}
index b2a5ae3..c28a0f3 100644 (file)
@@ -152,8 +152,8 @@ TString Setup
    //
    
    // add RSN task
-   gROOT->LoadMacro(Form("%s/AddAnalysisTaskRsnMini.C", macroPath));
-   if (!AddAnalysisTaskRsnMini(isMC, isPP, macroPath, nmix)) return "";
+   gROOT->LoadMacro(Form("%s/AddAnalysisTaskRsnMiniTest.C", macroPath));
+   if (!AddAnalysisTaskRsnMiniTest(isMC, isPP, macroPath, nmix)) return "";
    
    ::Info("AnalysisSetup", "Setup successful");
    return out;
diff --git a/PWG2/RESONANCES/macros/mini/ConfigKStarSimple.C b/PWG2/RESONANCES/macros/mini/ConfigKStarSimple.C
new file mode 100644 (file)
index 0000000..4d62cb5
--- /dev/null
@@ -0,0 +1,89 @@
+//
+// *** Configuration script for phi->KK analysis with 2010 runs ***
+// 
+// A configuration script for RSN package needs to define the followings:
+//
+// (1) decay tree of each resonance to be studied, which is needed to select
+//     true pairs and to assign the right mass to all candidate daughters
+// (2) cuts at all levels: single daughters, tracks, events
+// (3) output objects: histograms or trees
+//
+Bool_t ConfigKStarSimple
+(  
+   AliRsnMiniAnalysisTask *task, 
+   Bool_t                  isMC, 
+   
+   const char             *name,
+   const char             *outType,
+   const char             *computationType,
+   char                    charge1,
+   char                    charge2,
+   Bool_t                  useIM,
+   
+   AliRsnCutSet           *cutsPair
+)
+{
+   // 
+   // -- Define track cuts -------------------------------------------------------------------------
+   //
+   
+   // integrated pion cut
+   AliRsnCutDaughterKStar2010PP *cutPi = new AliRsnCutDaughterKStar2010PP("cutPionForKStar", AliPID::kPion);
+   // cut set
+   AliRsnCutSet *cutSetPi = new AliRsnCutSet("setPionForKStar", AliRsnTarget::kDaughter);
+   cutSetPi->AddCut(cutPi);
+   cutSetPi->SetCutScheme(cutPi->GetName());
+   // add to task
+   Int_t iCutPi = task->AddTrackCuts(cutSetPi);
+   
+   // integrated kaon cut
+   AliRsnCutDaughterKStar2010PP *cutK = new AliRsnCutDaughterKStar2010PP("cutKaonForKStar", AliPID::kKaon);
+   // cut set
+   AliRsnCutSet *cutSetK = new AliRsnCutSet("setKaonForKStar", AliRsnTarget::kDaughter);
+   cutSetK->AddCut(cutK);
+   cutSetK->SetCutScheme(cutK->GetName());
+   // add to task
+   Int_t iCutK = task->AddTrackCuts(cutSetK);
+   
+   //
+   // -- Values ------------------------------------------------------------------------------------
+   //
+   
+   /* invariant mass   */ Int_t imID   = task->CreateValue(AliRsnMiniValue::kInvMass, kFALSE);
+   /* IM resolution    */ Int_t resID  = task->CreateValue(AliRsnMiniValue::kInvMassRes, kTRUE);
+   /* transv. momentum */ Int_t ptID   = task->CreateValue(AliRsnMiniValue::kPt, kFALSE);
+   /* centrality       */ Int_t centID = task->CreateValue(AliRsnMiniValue::kMult, kFALSE);
+   
+   //
+   // -- Create all needed outputs -----------------------------------------------------------------
+   //
+   
+   // create output
+   AliRsnMiniOutput *out = task->CreateOutput(name, outType, computationType);
+   
+   // selection settings
+   out->SetCutID(0, iCutK);
+   out->SetCutID(1, iCutPi);
+   
+   // daughter species
+   out->SetDaughter(0, AliRsnDaughter::kKaon);
+   out->SetDaughter(1, AliRsnDaughter::kPion);
+   out->SetCharge(0, charge1);
+   out->SetCharge(1, charge2);
+   
+   // resonance properties
+   out->SetMotherPDG(313);
+   out->SetMotherMass(0.896);
+   
+   // pair cuts
+   out->SetPairCuts(cutsPair);
+   
+   // axis X: invmass (or resolution)
+   if (useIM) 
+      out->AddAxis(imID, 90, 0.6, 1.5);
+   else
+      out->AddAxis(resID, 200, -0.02, 0.02);
+      
+   // axis Y: transverse momentum
+   out->AddAxis(ptID, 100, 0.0, 10.0);
+}
diff --git a/PWG2/RESONANCES/macros/mini/ConfigPhiSimple.C b/PWG2/RESONANCES/macros/mini/ConfigPhiSimple.C
new file mode 100644 (file)
index 0000000..bdee928
--- /dev/null
@@ -0,0 +1,82 @@
+//
+// *** Configuration script for phi->KK analysis with 2010 runs ***
+// 
+// A configuration script for RSN package needs to define the followings:
+//
+// (1) decay tree of each resonance to be studied, which is needed to select
+//     true pairs and to assign the right mass to all candidate daughters
+// (2) cuts at all levels: single daughters, tracks, events
+// (3) output objects: histograms or trees
+//
+Bool_t ConfigPhiSimple
+(  
+   AliRsnMiniAnalysisTask *task, 
+   Bool_t                  isMC, 
+   
+   const char             *name,
+   const char             *outType,
+   const char             *computationType,
+   char                    charge1,
+   char                    charge2,
+   Bool_t                  useIM,
+   
+   AliRsnCutSet           *cutsPair
+)
+{
+   // 
+   // -- Define track cuts -------------------------------------------------------------------------
+   //
+   
+   // integrated kaon cut
+   AliRsnCutKaonForPhi2010PP *cutStd = new AliRsnCutKaonForPhi2010PP("cutStdPP");
+   // cut set
+   AliRsnCutSet *cutSetStd = new AliRsnCutSet("kaonForPhi", AliRsnTarget::kDaughter);
+   cutSetStd->AddCut(cutStd);
+   cutSetStd->SetCutScheme(cutStd->GetName());
+   // add to task
+   Int_t icut = task->AddTrackCuts(cutSetStd);
+   
+   //
+   // -- Values ------------------------------------------------------------------------------------
+   //
+   
+   /* invariant mass   */ Int_t imID   = task->CreateValue(AliRsnMiniValue::kInvMass, kFALSE);
+   /* IM resolution    */ Int_t resID  = task->CreateValue(AliRsnMiniValue::kInvMassRes, kTRUE);
+   /* transv. momentum */ Int_t ptID   = task->CreateValue(AliRsnMiniValue::kPt, kFALSE);
+   /* centrality       */ Int_t centID = task->CreateValue(AliRsnMiniValue::kMult, kFALSE);
+   
+   //
+   // -- Create all needed outputs -----------------------------------------------------------------
+   //
+   
+   // create output object                          "HIST"   "PAIR" or "MIX" or "TRUE"
+   AliRsnMiniOutput *out = task->CreateOutput(name, outType, computationType);
+      
+   // cut IDs
+   out->SetCutID(0, icut);
+   out->SetCutID(1, icut);
+   
+   // daughter species   
+   out->SetDaughter(0, AliRsnDaughter::kKaon);
+   out->SetDaughter(1, AliRsnDaughter::kKaon);
+   out->SetCharge(0, charge1);
+   out->SetCharge(1, charge2);
+   
+   // resonance properties
+   out->SetMotherPDG(333);
+   out->SetMotherMass(1.019455);
+   
+   // pair cuts
+   out->SetPairCuts(cutsPair);
+   
+   // axis X: invmass (or resolution)
+   if (useIM) 
+      out->AddAxis(imID, 500, 0.9,  1.4);
+   else
+      out->AddAxis(resID, 200, -0.02, 0.02);
+   
+   // axis Y: transverse momentum
+   out->AddAxis(ptID, 100, 0.0, 10.0);
+   
+   return kTRUE;
+}
diff --git a/PWG2/RESONANCES/macros/mini/ConfigSigmaStar.C b/PWG2/RESONANCES/macros/mini/ConfigSigmaStar.C
new file mode 100644 (file)
index 0000000..b08e2db
--- /dev/null
@@ -0,0 +1,119 @@
+//
+// *** Configuration script for phi->KK analysis with 2010 runs ***
+// 
+// A configuration script for RSN package needs to define the followings:
+//
+// (1) decay tree of each resonance to be studied, which is needed to select
+//     true pairs and to assign the right mass to all candidate daughters
+// (2) cuts at all levels: single daughters, tracks, events
+// (3) output objects: histograms or trees
+//
+Bool_t ConfigSigmaStar
+(  
+   AliRsnMiniAnalysisTask *task, 
+   Bool_t                  isMC, 
+   const char             *suffix,
+   AliRsnCutSet           *cutsPair
+)
+{
+   // manage suffix
+   if (strlen(suffix) > 0) suffix = Form("_%s", suffix);
+   
+   // 
+   // -- Define track cuts -------------------------------------------------------------------------
+   //
+   
+   // integrated pion cut
+   AliRsnCutDaughterKStar2010PP *cutPi = new AliRsnCutDaughterKStar2010PP("cutPionForKStar", AliPID::kPion);
+   // cut set
+   AliRsnCutSet *cutSetPi = new AliRsnCutSet("setPionForKStar", AliRsnTarget::kDaughter);
+   cutSetPi->AddCut(cutPi);
+   cutSetPi->SetCutScheme(cutPi->GetName());
+   // add to task
+   Int_t iCutPi = task->AddTrackCuts(cutSetPi);
+   
+   // quality cuts
+   AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("qualityDaughterLambda");
+   
+   esdTrackCuts->SetAcceptKinkDaughters(0); // 0 = kFalse
+   //esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
+   esdTrackCuts->SetMaxChi2PerClusterTPC(4);
+   esdTrackCuts->SetMinNClustersTPC(70);
+   esdTrackCuts->SetRequireTPCRefit();
+   
+   // cut lambda
+   AliRsnCutV0 *cutLambda = new AliRsnCutV0("cutLambda", kLambda0);
+   cutLambda->SetESDtrackCuts(esdTrackCuts);
+   // cut set
+   AliRsnCutSet *cutSetLambda = new AliRsnCutSet("setLambda", AliRsnTarget::kDaughter);
+   cutSetLambda->AddCut(cutLambda);
+   cutSetLambda->SetCutScheme(cutLambda->GetName());
+   // add to task
+   Int_t iCutLambda = task->AddTrackCuts(cutSetLambda);
+   
+   // cut anti-AntiLambda
+   AliRsnCutV0 *cutAntiLambda = new AliRsnCutV0("cutAntiLambda", kLambda0Bar);
+   cutAntiLambda->SetESDtrackCuts(esdTrackCuts);
+   // cut set
+   AliRsnCutSet *cutSetAntiLambda = new AliRsnCutSet("setAntiLambda", AliRsnTarget::kDaughter);
+   cutSetAntiLambda->AddCut(cutAntiLambda);
+   cutSetAntiLambda->SetCutScheme(cutAntiLambda->GetName());
+   // add to task
+   Int_t iCutAntiLambda = task->AddTrackCuts(cutSetAntiLambda);
+   
+   //
+   // -- Values ------------------------------------------------------------------------------------
+   //
+   
+   /* invariant mass   */ Int_t imID   = task->CreateValue(AliRsnMiniValue::kInvMass, kFALSE);
+   /* IM resolution    */ Int_t resID  = task->CreateValue(AliRsnMiniValue::kInvMassRes, kTRUE);
+   /* transv. momentum */ Int_t ptID   = task->CreateValue(AliRsnMiniValue::kPt, kFALSE);
+   /* centrality       */ Int_t centID = task->CreateValue(AliRsnMiniValue::kMult, kFALSE);
+   
+   //
+   // -- Create all needed outputs -----------------------------------------------------------------
+   //
+   
+   // use an array for more compact writing, which are different on mixing and charges
+   // [0] = unlike
+   // [1] = mixing
+   // [2] = like ++
+   // [3] = like --
+   Bool_t   use     [12] = { 1         ,  1         ,  1             ,  1             ,  1         ,  1         ,  1             ,  1             ,  1         ,  1         ,  1             ,  1             };
+   Bool_t   useIM   [12] = { 1         ,  1         ,  1             ,  1             ,  1         ,  1         ,  1             ,  1             ,  1         ,  1         ,  1             ,  1             };
+   TString  name    [12] = {"SigmaP"   , "SigmaM"   , "ASigmaP"      , "ASigmaM"      , "SigmaPmix", "SigmaMmix", "ASigmaPmix"   , "ASigmaMmix"   , "SigmaPt"  , "SigmaMt"  , "ASigmaPt"     , "ASigmaMt"     };
+   TString  comp    [12] = {"PAIR"     , "PAIR"     , "PAIR"         , "PAIR"         , "MIX"      , "MIX"      , "MIX"          , "MIX"          , "TRUE"     , "TRUE"     , "TRUE"         , "TRUE"         };
+   TString  output  [12] = {"HIST"     , "HIST"     , "HIST"         , "HIST"         , "HIST"     , "HIST"     , "HIST"         , "HIST"         , "HIST"     , "HIST"     , "HIST"         , "HIST"         };
+   Char_t   charge1 [12] = {'0'        , '0'        , '0'            , '0'            , '0'        , '0'        , '0'            , '0'            , '0'        , '0'        , '0'            , '0'            };
+   Char_t   charge2 [12] = {'+'        , '-'        , '-'            , '+'            , '+'        , '-'        , '-'            , '+'            , '+'        , '-'        , '-'            , '+'            };
+   Int_t    cutID1  [12] = { iCutLambda,  iCutLambda,  iCutAntiLambda,  iCutAntiLambda,  iCutLambda,  iCutLambda,  iCutAntiLambda,  iCutAntiLambda,  iCutLambda,  iCutLambda,  iCutAntiLambda,  iCutAntiLambda};
+   Int_t    cutID2  [12] = { iCutPi    ,  iCutPi    ,  iCutPi        ,  iCutPi        ,  iCutPi    ,  iCutPi    ,  iCutPi        ,  iCutPi        ,  iCutPi    ,  iCutPi    ,  iCutPi        ,  iCutPi        };
+   Int_t    ipdg    [12] = { 3224      ,  3114      , -3224          , -3114          ,  3224      ,  3114      , -3224          , -3114          ,  3224      ,  3114      , -3224          , -3114          };
+   Double_t mass    [12] = { 1382.3    ,  1387.4    ,  1382.3        ,  1387.4        ,  1382.3    ,  1387.4    ,  1382.3        ,  1387.4        ,  1382.3    ,  1387.4    ,  1382.3        ,  1387.4        };
+   
+   for (Int_t i = 0; i < 12; i++) {
+      if (!use[i]) continue;
+      // create output
+      AliRsnMiniOutput *out = task->CreateOutput(Form("sigmastar_%s%s", name[i].Data(), suffix), output[i].Data(), comp[i].Data());
+      // selection settings
+      out->SetCutID(0, cutID1[i]);
+      out->SetCutID(1, cutID2[i]);
+      out->SetDaughter(0, AliRsnDaughter::kLambda);
+      out->SetDaughter(1, AliRsnDaughter::kPion);
+      out->SetCharge(0, charge1[i]);
+      out->SetCharge(1, charge2[i]);
+      out->SetMotherPDG(ipdg[i]);
+      out->SetMotherMass(mass[i]);
+      // pair cuts
+      out->SetPairCuts(cutsPair);
+      // axis X: invmass (or resolution)
+      if (useIM[i]) 
+         out->AddAxis(imID, 800, 1.2, 2.0);
+      else
+         out->AddAxis(resID, 200, -0.02, 0.02);
+      // axis Y: transverse momentum
+      out->AddAxis(ptID, 100, 0.0, 10.0);
+   }
+   
+   return kTRUE;
+}
diff --git a/PWG2/RESONANCES/macros/mini/ConfigSigmaStarMC.C b/PWG2/RESONANCES/macros/mini/ConfigSigmaStarMC.C
new file mode 100644 (file)
index 0000000..25b212c
--- /dev/null
@@ -0,0 +1,72 @@
+//
+// *** Configuration script for phi->KK analysis with 2010 runs ***
+// 
+// A configuration script for RSN package needs to define the followings:
+//
+// (1) decay tree of each resonance to be studied, which is needed to select
+//     true pairs and to assign the right mass to all candidate daughters
+// (2) cuts at all levels: single daughters, tracks, events
+// (3) output objects: histograms or trees
+//
+Bool_t ConfigSigmaStarMC
+(
+   AliRsnMiniAnalysisTask *task, 
+   Bool_t                  isPP, 
+   const char             *suffix,
+   AliRsnCutSet           *cutsPair
+)
+{
+   // manage suffix
+   if (strlen(suffix) > 0) suffix = Form("_%s", suffix);
+   
+   // 
+   // -- Define track cuts -------------------------------------------------------------------------
+   //
+   
+   /*** EMPTY FOR TRUE PAIRS COMPUTATION ***/
+   
+   //
+   // -- Values ------------------------------------------------------------------------------------
+   //
+   
+   /* invariant mass   */ Int_t imID   = task->CreateValue(AliRsnMiniValue::kInvMass, kFALSE);
+   /* transv. momentum */ Int_t ptID   = task->CreateValue(AliRsnMiniValue::kPt, kFALSE);
+   /* centrality       */ Int_t centID = task->CreateValue(AliRsnMiniValue::kMult, kFALSE);
+
+   //
+   // -- Create all needed outputs -----------------------------------------------------------------
+   //
+   
+   TString mode = "HIST";
+   if (!isPP) mode = "SPARSE";
+   
+   // create output
+   AliRsnMiniOutput *out = task->CreateOutput(Form("sigmastarP_TrueMC%s", suffix), mode.Data(), "MOTHER");
+   // selection settings
+   out->SetDaughter(0, AliRsnDaughter::kLambda);
+   out->SetDaughter(1, AliRsnDaughter::kPion);
+   out->SetMotherPDG(3224);
+   out->SetMotherMass(1382.3);
+   // pair cuts
+   out->SetPairCuts(cutsPair);
+   // binnings
+   out->AddAxis(imID, 800, 1.2, 2.0);
+   out->AddAxis(ptID, 100, 0.0, 10.0);
+   if (!isPP) out->AddAxis(centID, 100, 0.0, 100.0);
+   
+   // create output
+   AliRsnMiniOutput *out = task->CreateOutput(Form("sigmastarM_TrueMC%s", suffix), mode.Data(), "MOTHER");
+   // selection settings
+   out->SetDaughter(0, AliRsnDaughter::kLambda);
+   out->SetDaughter(1, AliRsnDaughter::kPion);
+   out->SetMotherPDG(3114);
+   out->SetMotherMass(1387.4);
+   // pair cuts
+   out->SetPairCuts(cutsPair);
+   // binnings
+   out->AddAxis(imID, 800, 1.2, 2.0);
+   out->AddAxis(ptID, 100, 0.0, 10.0);
+   if (!isPP) out->AddAxis(centID, 100, 0.0, 100.0);
+   
+   return kTRUE;
+}
index 2d6e6e5..373de75 100644 (file)
@@ -1,6 +1,6 @@
 void runLocal
 (
-   Int_t       nReadFiles      =  0,
+   Int_t       nReadFiles      =  5,
    Int_t       nSkipFiles      =  0,
    Int_t       nmix            =  2,
    //const char *inputSource     = "file-collections/AOD048_LHC11a10b.txt",
index 01a3bd8..a4ebeb5 100644 (file)
@@ -1,6 +1,6 @@
 void runPluginProof
 (
-   const char *runMode     = "full",
+   const char *runMode     = "test",
    
    Int_t       nmix        = 50,
    
@@ -20,7 +20,7 @@ void runPluginProof
    
    const char *outName     = "proof.root",
    const char *macroPath   = ".",
-   const char *testFile    = "",
+   const char *testFile    = "file-collections/ESD_LHC10d1.txt",
    const char *addLibs     = "libEventMixing.so PWG2resonances.par",
    const char *addPars     = ""
 )