Adding new cut monitors
authorakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 12 Jun 2008 10:08:23 +0000 (10:08 +0000)
committerakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 12 Jun 2008 10:08:23 +0000 (10:08 +0000)
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.cxx [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.h [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.cxx [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.h [new file with mode: 0644]
PWG2/PWG2femtoscopyLinkDef.h
PWG2/libPWG2femtoscopy.pkg

diff --git a/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.cxx b/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.cxx
new file mode 100644 (file)
index 0000000..df3b752
--- /dev/null
@@ -0,0 +1,83 @@
+////////////////////////////////////////////////////////////////////////////////
+//                                                                            //
+// AliFemtoCutMonitorEventMult - the cut monitor for particles to study    //
+// the difference between reconstructed and true momentum                     //
+//                                                                            //
+////////////////////////////////////////////////////////////////////////////////
+#include "AliFemtoCutMonitorEventMult.h"
+#include "AliFemtoModelHiddenInfo.h"
+#include "AliFemtoEvent.h"
+#include <TH1D.h>
+#include <TH2D.h>
+#include <TList.h>
+
+AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult():
+  fEvMult(0)
+{
+  // Default constructor
+  fEvMult = new TH1D("EvMult", "Event Multiplicity", 5001, -0.5, 5000.5);
+}
+
+AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult(const char *aName):
+  AliFemtoCutMonitor(),
+  fEvMult(0)
+{
+  // Normal constructor
+  char name[200];
+  snprintf(name, 200, "EvMult%s", aName);
+  fEvMult = new TH1D(name, "Event Multiplicity", 5001, -0.5, 5000.5);
+}
+
+AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult(const AliFemtoCutMonitorEventMult &aCut):
+  AliFemtoCutMonitor(),
+  fEvMult(0)
+{
+  // copy constructor
+  if (fEvMult) delete fEvMult;
+  fEvMult = new TH1D(*aCut.fEvMult);
+}
+
+AliFemtoCutMonitorEventMult::~AliFemtoCutMonitorEventMult()
+{
+  // Destructor
+  delete fEvMult;
+}
+
+AliFemtoCutMonitorEventMult& AliFemtoCutMonitorEventMult::operator=(const AliFemtoCutMonitorEventMult& aCut)
+{
+  // assignment operator
+  if (this == &aCut) 
+    return *this;
+
+  if (fEvMult) delete fEvMult;
+  fEvMult = new TH1D(*aCut.fEvMult);
+  
+  return *this;
+}
+
+AliFemtoString AliFemtoCutMonitorEventMult::Report(){ 
+  // Prepare report from the execution
+  string stemp = "*** AliFemtoCutMonitorEventMult report"; 
+  AliFemtoString returnThis = stemp;
+  return returnThis; 
+}
+
+void AliFemtoCutMonitorEventMult::Fill(const AliFemtoEvent* aEvent)
+{
+  // Fill in the monitor histograms with the values from the current track
+  fEvMult->Fill(aEvent->NumberOfTracks());
+}
+
+void AliFemtoCutMonitorEventMult::Write()
+{
+  // Write out the relevant histograms
+  fEvMult->Write();
+}
+
+TList *AliFemtoCutMonitorEventMult::GetOutputList()
+{
+  TList *tOutputList = new TList();
+  tOutputList->Add(fEvMult);
+
+  return tOutputList;
+}
diff --git a/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.h b/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.h
new file mode 100644 (file)
index 0000000..b961ebe
--- /dev/null
@@ -0,0 +1,50 @@
+////////////////////////////////////////////////////////////////////////////////
+///                                                                          ///
+/// AliFemtoCutMonitorEventMult - the cut monitor for events to study        ///
+/// the multiplicity distribution of events                                  ///
+///                                                                          ///
+////////////////////////////////////////////////////////////////////////////////
+#ifndef AliFemtoCutMonitorEventMult_hh
+#define AliFemtoCutMonitorEventMult_hh
+
+class AliFemtoEvent;
+class AliFemtoTrack;
+class AliFemtoV0;
+class AliFemtoKink;
+class AliFemtoPair; 
+class TH1D;
+class TH2D;
+class TList;
+#include "AliFemtoString.h"
+#include "AliFemtoParticleCollection.h"
+#include "AliFemtoCutMonitor.h"
+
+class AliFemtoCutMonitorEventMult : public AliFemtoCutMonitor{
+  
+public:
+  AliFemtoCutMonitorEventMult();
+  AliFemtoCutMonitorEventMult(const char *aName);
+  AliFemtoCutMonitorEventMult(const AliFemtoCutMonitorEventMult &aCut);
+  virtual ~AliFemtoCutMonitorEventMult();
+
+  AliFemtoCutMonitorEventMult& operator=(const AliFemtoCutMonitorEventMult& aCut);
+
+  virtual AliFemtoString Report();
+  virtual void Fill(const AliFemtoEvent* aEvent);
+  virtual void Fill(const AliFemtoTrack* aTrack) {AliFemtoCutMonitor::Fill(aTrack);}
+  virtual void Fill(const AliFemtoV0* aV0) {AliFemtoCutMonitor::Fill(aV0);}
+  virtual void Fill(const AliFemtoKink* aKink) {AliFemtoCutMonitor::Fill(aKink);}
+  virtual void Fill(const AliFemtoPair* aPair) {AliFemtoCutMonitor::Fill(aPair);}
+  virtual void Fill(const AliFemtoParticleCollection* aCollection) {AliFemtoCutMonitor::Fill(aCollection);}
+  virtual void Fill(const AliFemtoEvent* aEvent,const AliFemtoParticleCollection* aCollection)
+  {AliFemtoCutMonitor::Fill(aEvent, aCollection);}
+
+  void Write();
+
+  virtual TList *GetOutputList();
+
+private:
+  TH1D *fEvMult;    // Multiplicity distribution
+};
+
+#endif
diff --git a/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.cxx b/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.cxx
new file mode 100644 (file)
index 0000000..d600a78
--- /dev/null
@@ -0,0 +1,103 @@
+////////////////////////////////////////////////////////////////////////////////
+//                                                                            //
+// AliFemtoCutMonitorParticleVertPos - the cut monitor for particles to study    //
+// the difference between reconstructed and true momentum                     //
+//                                                                            //
+////////////////////////////////////////////////////////////////////////////////
+#include "AliFemtoCutMonitorParticleVertPos.h"
+#include "AliFemtoModelHiddenInfo.h"
+#include <TH1D.h>
+#include <TH2D.h>
+#include <TList.h>
+#include <TMath.h>
+
+AliFemtoCutMonitorParticleVertPos::AliFemtoCutMonitorParticleVertPos():
+  fVertPos(0),
+  fEtaZ(0)
+{
+  // Default constructor
+  fVertPos = new TH2D("VertPos", "Vertex position", 200, -20.0, 20.0, 200, -20.0, 20.0);
+  fEtaZ    = new TH2D("EtaZPos", "Z vs. Eta", 200, -100.0, 100.0, 100, -1.5, 1.5);
+}
+
+AliFemtoCutMonitorParticleVertPos::AliFemtoCutMonitorParticleVertPos(const char *aName):
+  AliFemtoCutMonitor(),
+  fVertPos(0),
+  fEtaZ(0)
+{
+  // Normal constructor
+  char name[200];
+  snprintf(name, 200, "VertPos%s", aName);
+  fVertPos = new TH2D(name, "Rapdity vs Pt", 200, -20.0, 20.0, 200, -20.0, 20.0);
+  snprintf(name, 200, "EtaZPos%s", aName);
+  fEtaZ    = new TH2D(name, "Z vs. Eta", 200, -100.0, 100.0, 100, -1.5, 1.5);
+}
+
+AliFemtoCutMonitorParticleVertPos::AliFemtoCutMonitorParticleVertPos(const AliFemtoCutMonitorParticleVertPos &aCut):
+  AliFemtoCutMonitor(),
+  fVertPos(0),
+  fEtaZ(0)
+
+{
+  // copy constructor
+  if (fVertPos) delete fVertPos;
+  fVertPos = new TH2D(*aCut.fVertPos);
+  if (fEtaZ) delete fEtaZ;
+  fEtaZ = new TH2D(*aCut.fEtaZ);
+}
+
+AliFemtoCutMonitorParticleVertPos::~AliFemtoCutMonitorParticleVertPos()
+{
+  // Destructor
+  delete fVertPos;
+  delete fEtaZ;
+}
+
+AliFemtoCutMonitorParticleVertPos& AliFemtoCutMonitorParticleVertPos::operator=(const AliFemtoCutMonitorParticleVertPos& aCut)
+{
+  // assignment operator
+  if (this == &aCut) 
+    return *this;
+
+  if (fVertPos) delete fVertPos;
+  fVertPos = new TH2D(*aCut.fVertPos);
+  if (fEtaZ) delete fEtaZ;
+  fEtaZ = new TH2D(*aCut.fEtaZ);
+  
+  return *this;
+}
+
+AliFemtoString AliFemtoCutMonitorParticleVertPos::Report(){ 
+  // Prepare report from the execution
+  string stemp = "*** AliFemtoCutMonitorParticleVertPos report"; 
+  AliFemtoString returnThis = stemp;
+  return returnThis; 
+}
+
+void AliFemtoCutMonitorParticleVertPos::Fill(const AliFemtoTrack* aTrack)
+{
+  // Fill in the monitor histograms with the values from the current track
+  AliFemtoModelHiddenInfo *hinfo = dynamic_cast<AliFemtoModelHiddenInfo *>(aTrack->GetHiddenInfo());
+  if (hinfo) {
+    float tEta = -TMath::Log(TMath::Tan(hinfo->GetTrueMomentum()->theta()/2.0));
+
+    fVertPos->Fill(hinfo->GetEmissionPoint()->x(), hinfo->GetEmissionPoint()->y());
+    fEtaZ->Fill(hinfo->GetEmissionPoint()->z(), tEta);
+  }
+}
+
+void AliFemtoCutMonitorParticleVertPos::Write()
+{
+  // Write out the relevant histograms
+  fVertPos->Write();
+  fEtaZ->Write();
+}
+
+TList *AliFemtoCutMonitorParticleVertPos::GetOutputList()
+{
+  TList *tOutputList = new TList();
+  tOutputList->Add(fVertPos);
+  tOutputList->Add(fEtaZ);
+
+  return tOutputList;
+}
diff --git a/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.h b/PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.h
new file mode 100644 (file)
index 0000000..a185b29
--- /dev/null
@@ -0,0 +1,52 @@
+////////////////////////////////////////////////////////////////////////////////
+///                                                                          ///
+/// AliFemtoCutMonitorParticleVertPos - the cut monitor for particles to study  ///
+/// the difference between reconstructed and true momentum    ///
+///                                                                          ///
+////////////////////////////////////////////////////////////////////////////////
+#ifndef AliFemtoCutMonitorParticleVertPos_hh
+#define AliFemtoCutMonitorParticleVertPos_hh
+
+class AliFemtoEvent;
+class AliFemtoTrack;
+class AliFemtoV0;
+class AliFemtoKink;
+class AliFemtoPair; // Gael 12/04/02
+class TH1D;
+class TH2D;
+class TList;
+#include "AliFemtoString.h"
+#include "AliFemtoParticleCollection.h"
+#include "AliFemtoCutMonitor.h"
+
+class AliFemtoCutMonitorParticleVertPos : public AliFemtoCutMonitor{
+  
+public:
+  AliFemtoCutMonitorParticleVertPos();
+  AliFemtoCutMonitorParticleVertPos(const char *aName);
+  AliFemtoCutMonitorParticleVertPos(const AliFemtoCutMonitorParticleVertPos &aCut);
+  virtual ~AliFemtoCutMonitorParticleVertPos();
+
+  AliFemtoCutMonitorParticleVertPos& operator=(const AliFemtoCutMonitorParticleVertPos& aCut);
+
+  virtual AliFemtoString Report();
+  virtual void Fill(const AliFemtoEvent* aEvent) {AliFemtoCutMonitor::Fill(aEvent);}
+  virtual void Fill(const AliFemtoTrack* aTrack);
+  virtual void Fill(const AliFemtoV0* aV0) {AliFemtoCutMonitor::Fill(aV0);}
+  virtual void Fill(const AliFemtoKink* aKink) {AliFemtoCutMonitor::Fill(aKink);}
+  virtual void Fill(const AliFemtoPair* aPair) {AliFemtoCutMonitor::Fill(aPair);}
+  virtual void Fill(const AliFemtoParticleCollection* aCollection) {AliFemtoCutMonitor::Fill(aCollection);}
+  virtual void Fill(const AliFemtoEvent* aEvent,const AliFemtoParticleCollection* aCollection)
+  {AliFemtoCutMonitor::Fill(aEvent, aCollection);}
+
+
+  void Write();
+
+  virtual TList *GetOutputList();
+
+private:
+  TH2D *fVertPos;    // Vertex position x vs. y monitor
+  TH2D *fEtaZ;       // Vertex z position vs. eta monitor
+};
+
+#endif
index 727fc70..36d3c56 100644 (file)
@@ -46,6 +46,8 @@
 #pragma link C++ class AliFemtoModelCorrFctn+;
 #pragma link C++ class AliFemtoModelWeightGeneratorLednicky+;
 #pragma link C++ class AliFemtoCutMonitorParticleYPt+;
+#pragma link C++ class AliFemtoCutMonitorParticleVertPos+;
+#pragma link C++ class AliFemtoCutMonitorEventMult+;
 #pragma link C++ class AliFemtoEventReaderAOD+;
 #pragma link C++ class AliFemtoEventReaderAODChain+;
 #pragma link C++ class AliFemtoAODTrackCut+;
index 2f0bf1a..1880d87 100644 (file)
@@ -37,6 +37,8 @@ SRCS= FEMTOSCOPY/AliFemto/AliFemtoSimpleAnalysis.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoModelFreezeOutGenerator.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoModelWeightGeneratorLednicky.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleYPt.cxx \
+      FEMTOSCOPY/AliFemto/AliFemtoCutMonitorParticleVertPos.cxx \
+      FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoKTPairCut.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoCorrFctnNonIdDR.cxx \
       FEMTOSCOPY/AliFemto/AliFemtoCorrFctn3DSpherical.cxx \