Add Gamma conversions removal
authorakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 19 Apr 2010 08:21:57 +0000 (08:21 +0000)
committerakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 19 Apr 2010 08:21:57 +0000 (08:21 +0000)
PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.cxx [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.h [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.cxx [new file with mode: 0644]
PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.h [new file with mode: 0644]
PWG2/libPWG2femtoscopyUser.pkg

diff --git a/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.cxx b/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.cxx
new file mode 100644 (file)
index 0000000..611db02
--- /dev/null
@@ -0,0 +1,208 @@
+////////////////////////////////////////////////////////////////////////////////
+//                                                                            //
+// AliFemtoGammaMonitor - A correlation function that analyzes            //
+// two particle correlations with respect to the azimuthal angle (phi)        //
+// and pseudorapidity (eta) difference                                        //
+//                                                                            //
+// Authors: Adam Kisiel Adam.Kisiel@cern.ch                                   //
+//                                                                            //
+////////////////////////////////////////////////////////////////////////////////
+
+#include "AliFemtoCorrFctnGammaMonitor.h"
+//#include "AliFemtoHisto.hh"
+#include <cstdio>
+#include <TMath.h>
+
+#ifdef __ROOT__ 
+ClassImp(AliFemtoCorrFctnGammaMonitor)
+#endif
+
+//____________________________
+AliFemtoCorrFctnGammaMonitor::AliFemtoCorrFctnGammaMonitor(char* title, const int& aMinvBins=20, const int& aDThetaBins=20):
+  AliFemtoCorrFctn(),
+  fNumPMinvDTheta(0),
+  fDenPMinvDTheta(0),
+  fNumNMinvDTheta(0),
+  fDenNMinvDTheta(0)
+{
+  // set up numerator
+  char tTitNumD[100] = "NumPMinvTheta";
+  strcat(tTitNumD,title);
+  fNumPMinvDTheta = new TH2D(tTitNumD,title,aMinvBins,0.0,0.2,aDThetaBins,0.0,0.2);
+  // set up denominator
+  char tTitDenD[100] = "DenPMinvTheta";
+  strcat(tTitDenD,title);
+  fDenPMinvDTheta = new TH2D(tTitDenD,title,aMinvBins,0.0,0.2,aDThetaBins,0.0,0.2);
+
+  // set up numerator
+  char tTitNumR[100] = "NumNMinvTheta";
+  strcat(tTitNumR,title);
+  fNumNMinvDTheta = new TH2D(tTitNumR,title,aMinvBins,0.0,0.2,aDThetaBins,0.0,0.2);
+  // set up denominator
+  char tTitDenR[100] = "DenNMinvTheta";
+  strcat(tTitDenR,title);
+  fDenNMinvDTheta = new TH2D(tTitDenR,title,aMinvBins,0.0,0.2,aDThetaBins,0.0,0.2);
+
+  // to enable error bar calculation...
+  fNumPMinvDTheta->Sumw2();
+  fDenPMinvDTheta->Sumw2();
+  fNumNMinvDTheta->Sumw2();
+  fDenNMinvDTheta->Sumw2();
+}
+
+//____________________________
+AliFemtoCorrFctnGammaMonitor::AliFemtoCorrFctnGammaMonitor(const AliFemtoCorrFctnGammaMonitor& aCorrFctn) :
+  AliFemtoCorrFctn(),
+  fNumPMinvDTheta(0),
+  fDenPMinvDTheta(0),
+  fNumNMinvDTheta(0),
+  fDenNMinvDTheta(0)
+{
+  // copy constructor
+  if (aCorrFctn.fNumPMinvDTheta)
+    fNumPMinvDTheta = new TH2D(*aCorrFctn.fNumPMinvDTheta);
+  else
+    fNumPMinvDTheta = 0;
+  if (aCorrFctn.fDenPMinvDTheta)
+    fDenPMinvDTheta = new TH2D(*aCorrFctn.fDenPMinvDTheta);
+  else
+    fDenPMinvDTheta = 0;
+
+  if (aCorrFctn.fNumNMinvDTheta)
+    fNumNMinvDTheta = new TH2D(*aCorrFctn.fNumNMinvDTheta);
+  else
+    fNumNMinvDTheta = 0;
+  if (aCorrFctn.fDenNMinvDTheta)
+    fDenNMinvDTheta = new TH2D(*aCorrFctn.fDenNMinvDTheta);
+  else
+    fDenNMinvDTheta = 0;
+
+}
+//____________________________
+AliFemtoCorrFctnGammaMonitor::~AliFemtoCorrFctnGammaMonitor(){
+  // destructor
+  delete fNumPMinvDTheta;
+  delete fDenPMinvDTheta;
+  delete fNumNMinvDTheta;
+  delete fDenNMinvDTheta;
+}
+//_________________________
+AliFemtoCorrFctnGammaMonitor& AliFemtoCorrFctnGammaMonitor::operator=(const AliFemtoCorrFctnGammaMonitor& aCorrFctn)
+{
+  // assignment operator
+  if (this == &aCorrFctn)
+    return *this;
+
+  if (aCorrFctn.fNumPMinvDTheta)
+    fNumPMinvDTheta = new TH2D(*aCorrFctn.fNumPMinvDTheta);
+  else
+    fNumPMinvDTheta = 0;
+  if (aCorrFctn.fDenPMinvDTheta)
+    fDenPMinvDTheta = new TH2D(*aCorrFctn.fDenPMinvDTheta);
+  else
+    fDenPMinvDTheta = 0;
+
+  if (aCorrFctn.fNumNMinvDTheta)
+    fNumNMinvDTheta = new TH2D(*aCorrFctn.fNumNMinvDTheta);
+  else
+    fNumNMinvDTheta = 0;
+  if (aCorrFctn.fDenNMinvDTheta)
+    fDenNMinvDTheta = new TH2D(*aCorrFctn.fDenNMinvDTheta);
+  else
+    fDenNMinvDTheta = 0;
+
+
+  return *this;
+}
+//_________________________
+void AliFemtoCorrFctnGammaMonitor::Finish(){
+  // here is where we should normalize, fit, etc...
+  // we should NOT Draw() the histos (as I had done it below),
+  // since we want to insulate ourselves from root at this level
+  // of the code.  Do it instead at root command line with browser.
+  //  mShareNumerator->Draw();
+  //mShareDenominator->Draw();
+  //mRatio->Draw();
+
+}
+
+//____________________________
+AliFemtoString AliFemtoCorrFctnGammaMonitor::Report(){
+  // create report
+  string stemp = "Gamma Monitor Function Report:\n";
+  char ctemp[100];
+  sprintf(ctemp,"Number of entries in numerator:\t%E\n",fNumPMinvDTheta->GetEntries());
+  stemp += ctemp;
+  sprintf(ctemp,"Number of entries in denominator:\t%E\n",fDenPMinvDTheta->GetEntries());
+  stemp += ctemp;
+  //  stemp += mCoulombWeight->Report();
+  AliFemtoString returnThis = stemp;
+  return returnThis;
+}
+//____________________________
+void AliFemtoCorrFctnGammaMonitor::AddRealPair( AliFemtoPair* pair){
+  // add real (effect) pair
+  double me = 0.000511;
+
+  double theta1 = pair->Track1()->Track()->P().Theta();
+  double theta2 = pair->Track2()->Track()->P().Theta();
+  double dtheta = TMath::Abs(theta1 - theta2);
+
+  double e1 = TMath::Sqrt(me*me + pair->Track1()->Track()->P().Mag2());
+  double e2 = TMath::Sqrt(me*me + pair->Track2()->Track()->P().Mag2());
+
+  double minv = 2*me*me + 2*(e1*e2 - 
+                            pair->Track1()->Track()->P().x()*pair->Track2()->Track()->P().x() -
+                            pair->Track1()->Track()->P().y()*pair->Track2()->Track()->P().y() -
+                            pair->Track1()->Track()->P().z()*pair->Track2()->Track()->P().z());
+
+  if (pair->KSide()>0.0) 
+    fNumPMinvDTheta->Fill(minv, dtheta);
+  else
+    fNumNMinvDTheta->Fill(minv, dtheta);
+}
+//____________________________
+void AliFemtoCorrFctnGammaMonitor::AddMixedPair( AliFemtoPair* pair){
+  // add mixed (background) pair
+  double me = 0.000511;
+
+  double theta1 = pair->Track1()->Track()->P().Theta();
+  double theta2 = pair->Track2()->Track()->P().Theta();
+  double dtheta = TMath::Abs(theta1 - theta2);
+
+  double e1 = TMath::Sqrt(me*me + pair->Track1()->Track()->P().Mag2());
+  double e2 = TMath::Sqrt(me*me + pair->Track2()->Track()->P().Mag2());
+
+  double minv = 2*me*me + 2*(e1*e2 - 
+                            pair->Track1()->Track()->P().x()*pair->Track2()->Track()->P().x() -
+                            pair->Track1()->Track()->P().y()*pair->Track2()->Track()->P().y() -
+                            pair->Track1()->Track()->P().z()*pair->Track2()->Track()->P().z());
+
+  if (pair->KSide()>0.0) 
+    fDenPMinvDTheta->Fill(minv, dtheta);
+  else
+    fDenNMinvDTheta->Fill(minv, dtheta);
+}
+
+
+void AliFemtoCorrFctnGammaMonitor::WriteHistos()
+{
+  // Write out result histograms
+  fNumPMinvDTheta->Write();
+  fDenPMinvDTheta->Write();
+  fNumNMinvDTheta->Write();
+  fDenNMinvDTheta->Write();
+}
+
+TList* AliFemtoCorrFctnGammaMonitor::GetOutputList()
+{
+  // Prepare the list of objects to be written to the output
+  TList *tOutputList = new TList();
+
+  tOutputList->Add(fNumPMinvDTheta);
+  tOutputList->Add(fDenPMinvDTheta);
+  tOutputList->Add(fNumNMinvDTheta);
+  tOutputList->Add(fDenNMinvDTheta);
+
+  return tOutputList;
+}
diff --git a/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.h b/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.h
new file mode 100644 (file)
index 0000000..dc1a62c
--- /dev/null
@@ -0,0 +1,49 @@
+////////////////////////////////////////////////////////////////////////////////
+//                                                                            //
+// AliFemtoCorrFctnGammaMonitor - A correlation function that analyzes            //
+// two particle correlations with respect to the azimuthal angle (phi)        //
+// and pseudorapidity (eta) difference                                        //
+//                                                                            //
+// Authors: Adam Kisiel Adam.Kisiel@cern.ch                                   //
+//                                                                            //
+////////////////////////////////////////////////////////////////////////////////
+
+#ifndef ALIFEMTOCORRFCTNGAMMAMONITOR_H
+#define ALIFEMTOCORRFCTNGAMMAMONITOR_H
+
+#include "TH1D.h"
+#include "TH2D.h"
+#include "AliFemtoCorrFctn.h"
+
+class AliFemtoCorrFctnGammaMonitor : public AliFemtoCorrFctn {
+public:
+  AliFemtoCorrFctnGammaMonitor(char* title, const int& aMinvBins, const int& aDThetaBins);
+  AliFemtoCorrFctnGammaMonitor(const AliFemtoCorrFctnGammaMonitor& aCorrFctn);
+  virtual ~AliFemtoCorrFctnGammaMonitor();
+
+  AliFemtoCorrFctnGammaMonitor& operator=(const AliFemtoCorrFctnGammaMonitor& aCorrFctn);
+
+  virtual AliFemtoString Report();
+  virtual void AddRealPair(AliFemtoPair* aPair);
+  virtual void AddMixedPair(AliFemtoPair* aPair);
+
+  virtual void Finish();
+
+  void WriteHistos();
+  virtual TList* GetOutputList();
+private:
+  
+  TH2D *fNumPMinvDTheta;        // Numerator Minv vs. DTheta Positive kSide
+  TH2D *fDenPMinvDTheta;        // Denominator Minv vs. DTheta Positive kSide
+
+  TH2D *fNumNMinvDTheta;        // Numerator Minv vs. DTheta Negative kSide
+  TH2D *fDenNMinvDTheta;        // Denominator Minv vs. DTheta Negative kSide
+
+#ifdef __ROOT__
+  ClassDef(AliFemtoCorrFctnGammaMonitor, 1)
+#endif
+};
+
+
+#endif
+
diff --git a/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.cxx b/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.cxx
new file mode 100644 (file)
index 0000000..42e7b0c
--- /dev/null
@@ -0,0 +1,121 @@
+/////////////////////////////////////////////////////////////////////////////
+//                                                                         //
+// AliFemtoPairCutAntiGamma - a pair cut which checks     //
+// for some pair qualities that attempt to identify slit/doubly            //
+// reconstructed tracks and also selects pairs based on their separation   //
+// at the entrance to the TPC                                              //
+//                                                                         //
+/////////////////////////////////////////////////////////////////////////////
+/***************************************************************************
+ *
+ * $Id: AliFemtoPairCutAntiGamma.cxx,v 1.1.2.1 2007/10/19 13:35:33 akisiel Exp $
+ *
+ * Author: Adam Kisiel, Ohio State, kisiel@mps.ohio-state.edu
+ ***************************************************************************
+ *
+ * Description: part of STAR HBT Framework: AliFemtoMaker package
+ *   a cut to remove "shared" and "split" pairs
+ *
+ ***************************************************************************
+ *
+ *
+ **************************************************************************/
+
+#include "AliFemtoPairCutAntiGamma.h"
+#include <string>
+#include <cstdio>
+#include <TMath.h>
+
+#ifdef __ROOT__
+ClassImp(AliFemtoPairCutAntiGamma)
+#endif
+
+//__________________
+AliFemtoPairCutAntiGamma::AliFemtoPairCutAntiGamma():
+  AliFemtoShareQualityPairCut(),
+  fMaxEEMinv(0.0),
+  fMaxDTheta(0.0)
+{
+}
+//__________________
+AliFemtoPairCutAntiGamma::AliFemtoPairCutAntiGamma(const AliFemtoPairCutAntiGamma& c) : 
+  AliFemtoShareQualityPairCut(c),
+  fMaxEEMinv(0.0),
+  fMaxDTheta(0.0)
+{ 
+  fMaxEEMinv = c.fMaxEEMinv;
+  fMaxDTheta = c.fMaxDTheta;
+}
+
+//__________________
+AliFemtoPairCutAntiGamma::~AliFemtoPairCutAntiGamma(){
+  /* no-op */
+}
+//__________________
+bool AliFemtoPairCutAntiGamma::Pass(const AliFemtoPair* pair){
+  // Accept pairs based on their TPC entrance separation and
+  // quality and sharity
+  bool temp = true;
+
+  double me = 0.000511;
+
+  if ((pair->Track1()->Track()->Charge() * pair->Track2()->Track()->Charge()) < 0.0) {
+    double theta1 = pair->Track1()->Track()->P().Theta();
+    double theta2 = pair->Track2()->Track()->P().Theta();
+    double dtheta = TMath::Abs(theta1 - theta2);
+    
+    double e1 = TMath::Sqrt(me*me + pair->Track1()->Track()->P().Mag2());
+    double e2 = TMath::Sqrt(me*me + pair->Track2()->Track()->P().Mag2());
+    
+    double minv = 2*me*me + 2*(e1*e2 - 
+                              pair->Track1()->Track()->P().x()*pair->Track2()->Track()->P().x() -
+                              pair->Track1()->Track()->P().y()*pair->Track2()->Track()->P().y() -
+                              pair->Track1()->Track()->P().z()*pair->Track2()->Track()->P().z());
+    
+    if ((minv < fMaxEEMinv) && (dtheta < fMaxDTheta)) temp = false;
+  }
+
+  if (temp) {
+    temp = AliFemtoShareQualityPairCut::Pass(pair);
+    if (temp) fNPairsPassed++;
+    else fNPairsFailed++;
+  }
+  else
+    fNPairsFailed++;
+
+  return temp;
+}
+//__________________
+AliFemtoString AliFemtoPairCutAntiGamma::Report(){
+  // Prepare a report from the execution
+  string stemp = "AliFemtoPairCutAntiGamma Pair Cut - remove pairs possibly coming from Gamma conversions\n";  
+  char ctemp[100];
+  stemp += ctemp;
+  sprintf(ctemp,"Number of pairs which passed:\t%ld  Number which failed:\t%ld\n",fNPairsPassed,fNPairsFailed);
+  stemp += ctemp;
+  AliFemtoString returnThis = stemp;
+  return returnThis;}
+//__________________
+
+TList *AliFemtoPairCutAntiGamma::ListSettings()
+{
+  // return a list of settings in a writable form
+  TList *tListSetttings =  AliFemtoShareQualityPairCut::ListSettings();
+  char buf[200];
+  snprintf(buf, 200, "AliFemtoPairCutAntiGamma.maxeeminv=%f", fMaxEEMinv);
+  snprintf(buf, 200, "AliFemtoPairCutAntiGamma.maxdtheta=%f", fMaxDTheta);
+  tListSetttings->AddLast(new TObjString(buf));
+
+  return tListSetttings;
+}
+
+void AliFemtoPairCutAntiGamma::SetMaxEEMinv(Double_t maxeeminv)
+{
+  fMaxEEMinv = maxeeminv;
+}
+
+void AliFemtoPairCutAntiGamma::SetMaxThetaDiff(Double_t maxdtheta)
+{
+  fMaxDTheta = maxdtheta;
+}
diff --git a/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.h b/PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.h
new file mode 100644 (file)
index 0000000..547e916
--- /dev/null
@@ -0,0 +1,55 @@
+/////////////////////////////////////////////////////////////////////////////
+//                                                                         //
+// AliFemtoPairCutAntiGamma - a pair cut which checks     //
+// for some pair qualities that attempt to identify slit/doubly            //
+// reconstructed tracks and also selects pairs based on their separation   //
+// at the entrance to the TPC                                              //
+//                                                                         //
+/////////////////////////////////////////////////////////////////////////////
+/***************************************************************************
+ *
+ * $Id: AliFemtoPairCutAntiGamma.h,v 1.1.2.1 2007/10/19 13:35:33 akisiel Exp $
+ *
+ * Author: Adam Kisiel, Ohio State University, kisiel@mps.ohio-state.edu
+ ***************************************************************************
+ *
+ * Description: part of STAR HBT Framework: AliFemtoMaker package
+ *   a cut to remove "shared" and "split" pairs
+ *
+ ***************************************************************************
+ *
+ *
+ **************************************************************************/
+
+
+#ifndef ALIFEMTOPAIRCUTANTIGAMMA_H
+#define ALIFEMTOPAIRCUTANTIGAMMA_H
+
+#include "AliFemtoPairCut.h"
+#include "AliFemtoShareQualityPairCut.h"
+
+class AliFemtoPairCutAntiGamma : public AliFemtoShareQualityPairCut{
+public:
+  AliFemtoPairCutAntiGamma();
+  AliFemtoPairCutAntiGamma(const AliFemtoPairCutAntiGamma& c);
+  virtual ~AliFemtoPairCutAntiGamma();
+
+  virtual bool Pass(const AliFemtoPair* pair);
+  virtual AliFemtoString Report();
+  virtual TList *ListSettings();
+  virtual AliFemtoPairCut* Clone();
+  void SetMaxEEMinv(Double_t maxeeminv);
+  void SetMaxThetaDiff(Double_t maxdtheta);
+  
+ protected:
+  Double_t fMaxEEMinv; // Maximum allowed ee Minv
+  Double_t fMaxDTheta; // Maximum polar angle difference
+
+#ifdef __ROOT__
+  ClassDef(AliFemtoPairCutAntiGamma, 0)
+#endif
+};
+
+inline AliFemtoPairCut* AliFemtoPairCutAntiGamma::Clone() { AliFemtoPairCutAntiGamma* c = new AliFemtoPairCutAntiGamma(*this); return c;}
+
+#endif
index 573d1d3..7be4e7a 100644 (file)
@@ -1,7 +1,6 @@
 #-*- Mode: Makefile -*-
 
-SRCS= FEMTOSCOPY/AliFemtoUser/AliFemtoQPairCut.cxx \
-      FEMTOSCOPY/AliFemtoUser/AliFemtoShareQualityPairCut.cxx \
+SRCS= FEMTOSCOPY/AliFemtoUser/AliFemtoShareQualityPairCut.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoShareQualityKTPairCut.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoShareQualityTPCEntranceSepPairCut.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoESDTrackCut.cxx \
@@ -30,6 +29,8 @@ SRCS= FEMTOSCOPY/AliFemtoUser/AliFemtoQPairCut.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoCutMonitorParticlePtPDG.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnTPCNcls.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnDEtaDPhi.cxx \
+      FEMTOSCOPY/AliFemtoUser/AliFemtoCorrFctnGammaMonitor.cxx \
+      FEMTOSCOPY/AliFemtoUser/AliFemtoPairCutAntiGamma.cxx \
       FEMTOSCOPY/AliFemtoUser/AliFemtoCutMonitorParticleEtCorr.cxx
 
 #      FEMTOSCOPY/AliFemtoUser/AliFemtoPhiPairCut.cxx \