Adding two track QA macros
authorakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 15 May 2009 16:00:06 +0000 (16:00 +0000)
committerakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 15 May 2009 16:00:06 +0000 (16:00 +0000)
PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/ConfigFemtoAnalysis.C [new file with mode: 0644]
PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/FemtoQAPlots.C [new file with mode: 0644]

diff --git a/PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/ConfigFemtoAnalysis.C b/PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/ConfigFemtoAnalysis.C
new file mode 100644 (file)
index 0000000..17414ed
--- /dev/null
@@ -0,0 +1,454 @@
+/*********************************************************************
+ *                                                                   *
+ * ConfigFemtoAnalysis.C - configuration macro for the femtoscopic   *
+ * analysis, meant as a QA process for two-particle effects          *
+ *                                                                   *
+ * Author: Adam Kisiel (Adam.Kisiel@cern.ch)                         *
+ *                                                                   *
+ *********************************************************************/
+
+#include "AliFemtoManager.h"
+#include "AliFemtoEventReaderESDChain.h"
+#include "AliFemtoEventReaderESDChainKine.h"
+#include "AliFemtoSimpleAnalysis.h"
+#include "AliFemtoBasicEventCut.h"
+#include "AliFemtoESDTrackCut.h"
+#include "AliFemtoCorrFctn.h"
+#include "AliFemtoCutMonitorParticleYPt.h"
+#include "AliFemtoCutMonitorParticleVertPos.h"
+#include "AliFemtoCutMonitorParticleMomRes.h"
+#include "AliFemtoCutMonitorEventMult.h"
+#include "AliFemtoCutMonitorEventVertex.h"
+#include "AliFemtoShareQualityTPCEntranceSepPairCut.h"
+#include "AliFemtoQinvCorrFctn.h"
+#include "AliFemtoShareQualityCorrFctn.h"
+#include "AliFemtoTPCInnerCorrFctn.h"
+#include "AliFemtoVertexMultAnalysis.h"
+#include "AliFemtoCorrFctn3DSpherical.h"
+#include "AliFemtoChi2CorrFctn.h"
+#include "AliFemtoCorrFctnTPCNcls.h"
+#include "AliFemtoModelBPLCMSCorrFctn.h"
+#include "AliFemtoModelCorrFctn3DSpherical.h"
+#include "AliFemtoModelGausLCMSFreezeOutGenerator.h"
+#include "AliFemtoModelGausRinvFreezeOutGenerator.h"
+#include "AliFemtoModelManager.h"
+#include "AliFemtoModelWeightGeneratorBasic.h"
+#include "AliFemtoModelWeightGeneratorLednicky.h"
+#include "AliFemtoCorrFctnDirectYlm.h"
+#include "AliFemtoModelCorrFctnDirectYlm.h"
+#include "AliFemtoModelCorrFctnSource.h"
+#include "AliFemtoCutMonitorParticlePtPDG.h"
+#include "AliFemtoKTPairCut.h"
+
+//________________________________________________________________________
+AliFemtoManager* ConfigFemtoAnalysis() {
+
+  double PionMass = 0.13956995;
+  double KaonMass = 0.493677;
+  
+//   AliFemtoEventReaderESDChainKine* Reader=new AliFemtoEventReaderESDChainKine();
+//   Reader->SetConstrained(true);
+//   Reader->SetUseTPCOnly(false);
+
+  AliFemtoEventReaderESDChain* Reader=new AliFemtoEventReaderESDChain();
+  Reader->SetConstrained(true);
+  Reader->SetUseTPCOnly(false);
+
+  AliFemtoManager* Manager=new AliFemtoManager();
+  Manager->SetEventReader(Reader);
+       
+  // Switches for QA analyses
+  int runPiPlusStandard = 1;
+  int runPiMinusStandard = 1;
+  int runPositiveQA = 1;
+  int runNegativeQA = 1;
+  int runPositiveTPCQA = 1;
+  int runNegativeTPCQA = 1;
+
+  // *** First QA task - standard HBT analysis with all the cut on ***
+
+  // *** Begin pion-pion (positive) analysis ***
+  if (runPiPlusStandard) {
+    AliFemtoVertexMultAnalysis *anpipstd = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpipstd->SetNumEventsToMix(10);
+    anpipstd->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpipstd = new AliFemtoBasicEventCut();
+    mecpipstd->SetEventMult(2,100000);
+    mecpipstd->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpipstd = new AliFemtoESDTrackCut();
+    dtcpipstd->SetPidProbPion(0.2,1.001);
+    dtcpipstd->SetPidProbMuon(0.0,1.0);
+    dtcpipstd->SetPidProbKaon(0.0,1.0);
+    dtcpipstd->SetPidProbProton(0.0,1.0);
+    dtcpipstd->SetCharge(1.0);
+    dtcpipstd->SetPt(0.15,0.5);
+    dtcpipstd->SetMass(PionMass);
+    // Track quality cuts
+    dtcpipstd->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    //  dtcpipstd->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpipstd->SetminTPCncls(50);
+    dtcpipstd->SetRemoveKinks(kTRUE);
+    dtcpipstd->SetLabel(kFALSE);
+    dtcpipstd->SetMaxITSChiNdof(2.5);
+    dtcpipstd->SetMaxTPCChiNdof(3.0);
+    dtcpipstd->SetMaxImpactXY(3.0);
+    dtcpipstd->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpipstd = new AliFemtoCutMonitorParticleYPt("cutPasspipstd", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpipstd = new AliFemtoCutMonitorParticleYPt("cutFailpipstd", 0.13957);
+    dtcpipstd->AddCutMonitor(cutPassYPtpipstd, cutFailYPtpipstd);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpipstd = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpipstd->SetShareQualityMax(0.0);
+    sqpcpipstd->SetShareFractionMax(0.02);
+    sqpcpipstd->SetRemoveSameLabel(kFALSE);
+    sqpcpipstd->SetTPCEntranceSepMinimum(2.0);
+
+    anpipstd->SetEventCut(mecpipstd);
+    anpipstd->SetFirstParticleCut(dtcpipstd);
+    anpipstd->SetSecondParticleCut(dtcpipstd);
+    anpipstd->SetPairCut(sqpcpipstd);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpipstd= new AliFemtoShareQualityCorrFctn("sqqinvcfpipstd",40,0.0,0.4);
+    anpipstd->AddCorrFctn(csqqinvpipstd);
+
+    AliFemtoCorrFctnDirectYlm *cylmpipstd = new AliFemtoCorrFctnDirectYlm("cylmpipstd",3,60,0.0,0.3);
+    anpipstd->AddCorrFctn(cylmpipstd);
+    
+    AliFemtoQinvCorrFctn *cqinvpipstd = new AliFemtoQinvCorrFctn("qinvcfpipstd", 100,0.0,1.0);
+    anpipstd->AddCorrFctn(cqinvpipstd);
+
+    AliFemtoChi2CorrFctn *cchiqinvpipstd= new AliFemtoChi2CorrFctn("chicfpipstd",40,0.0,0.4);
+    anpipstd->AddCorrFctn(cchiqinvpipstd);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspipstd = new AliFemtoCorrFctnTPCNcls("cqtpcnclspipstd",40,0.0,0.4);
+    anpipstd->AddCorrFctn(cqtpcnclspipstd);
+
+    Manager->AddAnalysis(anpipstd);    
+  }
+  // *** End pion-pion (positive) analysis
+
+  // *** Begin pion-pion (negative) analysis ***
+  if (runPiMinusStandard) {
+    AliFemtoVertexMultAnalysis *anpimstd = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpimstd->SetNumEventsToMix(10);
+    anpimstd->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpimstd = new AliFemtoBasicEventCut();
+    mecpimstd->SetEventMult(2,100000);
+    mecpimstd->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpimstd = new AliFemtoESDTrackCut();
+    dtcpimstd->SetPidProbPion(0.2,1.001);
+    dtcpimstd->SetPidProbMuon(0.0,1.0);
+    dtcpimstd->SetPidProbKaon(0.0,1.0);
+    dtcpimstd->SetPidProbProton(0.0,1.0);
+    dtcpimstd->SetCharge(-1.0);
+    dtcpimstd->SetPt(0.15,0.5);
+    dtcpimstd->SetMass(PionMass);
+    // Track quality cuts
+    dtcpimstd->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    //  dtcpimstd->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpimstd->SetminTPCncls(50);
+    dtcpimstd->SetRemoveKinks(kTRUE);
+    dtcpimstd->SetLabel(kFALSE);
+    dtcpimstd->SetMaxITSChiNdof(2.5);
+    dtcpimstd->SetMaxTPCChiNdof(3.0);
+    dtcpimstd->SetMaxImpactXY(3.0);
+    dtcpimstd->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpimstd = new AliFemtoCutMonitorParticleYPt("cutPasspimstd", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpimstd = new AliFemtoCutMonitorParticleYPt("cutFailpimstd", 0.13957);
+    dtcpimstd->AddCutMonitor(cutPassYPtpimstd, cutFailYPtpimstd);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimstd = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpimstd->SetShareQualityMax(0.0);
+    sqpcpimstd->SetShareFractionMax(0.02);
+    sqpcpimstd->SetRemoveSameLabel(kFALSE);
+    sqpcpimstd->SetTPCEntranceSepMinimum(2.0);
+
+    anpimstd->SetEventCut(mecpimstd);
+    anpimstd->SetFirstParticleCut(dtcpimstd);
+    anpimstd->SetSecondParticleCut(dtcpimstd);
+    anpimstd->SetPairCut(sqpcpimstd);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpimstd= new AliFemtoShareQualityCorrFctn("sqqinvcfpimstd",40,0.0,0.4);
+    anpimstd->AddCorrFctn(csqqinvpimstd);
+
+    AliFemtoCorrFctnDirectYlm *cylmpimstd = new AliFemtoCorrFctnDirectYlm("cylmpimstd",3,60,0.0,0.3);
+    anpimstd->AddCorrFctn(cylmpimstd);
+    
+    AliFemtoQinvCorrFctn *cqinvpimstd = new AliFemtoQinvCorrFctn("qinvcfpimstd", 100,0.0,1.0);
+    anpimstd->AddCorrFctn(cqinvpimstd);
+
+    AliFemtoChi2CorrFctn *cchiqinvpimstd= new AliFemtoChi2CorrFctn("chicfpimstd",40,0.0,0.4);
+    anpimstd->AddCorrFctn(cchiqinvpimstd);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspimstd = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimstd",40,0.0,0.4);
+    anpimstd->AddCorrFctn(cqtpcnclspimstd);
+
+    Manager->AddAnalysis(anpimstd);    
+  }
+  // *** End pion-pion (negative) analysis
+
+  // *** Second QA task - HBT analysis with all pair cuts off ***
+  // *** Begin pion-pion (positive) analysis ***
+  if (runPositiveQA) {
+    AliFemtoVertexMultAnalysis *anpipnct = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpipnct->SetNumEventsToMix(10);
+    anpipnct->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpipnct = new AliFemtoBasicEventCut();
+    mecpipnct->SetEventMult(2,100000);
+    mecpipnct->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpipnct = new AliFemtoESDTrackCut();
+    dtcpipnct->SetPidProbPion(0.0,1.001);
+    dtcpipnct->SetPidProbMuon(0.0,1.0);
+    dtcpipnct->SetPidProbKaon(0.0,1.0);
+    dtcpipnct->SetPidProbProton(0.0,1.0);
+    dtcpipnct->SetCharge(1.0);
+    dtcpipnct->SetPt(0.1,1.0);
+    dtcpipnct->SetMass(PionMass);
+    // Track quality cuts
+    dtcpipnct->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    //  dtcpipnct->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpipnct->SetminTPCncls(50);
+    dtcpipnct->SetRemoveKinks(kTRUE);
+    dtcpipnct->SetLabel(kFALSE);
+    dtcpipnct->SetMaxITSChiNdof(6.0);
+    dtcpipnct->SetMaxTPCChiNdof(6.0);
+    dtcpipnct->SetMaxImpactXY(3.0);
+    dtcpipnct->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpipnct = new AliFemtoCutMonitorParticleYPt("cutPasspipnct", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpipnct = new AliFemtoCutMonitorParticleYPt("cutFailpipnct", 0.13957);
+    dtcpipnct->AddCutMonitor(cutPassYPtpipnct, cutFailYPtpipnct);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpipnct = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpipnct->SetShareQualityMax(1.0);
+    sqpcpipnct->SetShareFractionMax(1.0);
+    sqpcpipnct->SetRemoveSameLabel(kFALSE);
+    sqpcpipnct->SetTPCEntranceSepMinimum(0.0);
+
+    anpipnct->SetEventCut(mecpipnct);
+    anpipnct->SetFirstParticleCut(dtcpipnct);
+    anpipnct->SetSecondParticleCut(dtcpipnct);
+    anpipnct->SetPairCut(sqpcpipnct);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpipnct= new AliFemtoShareQualityCorrFctn("sqqinvcfpipnct",40,0.0,0.4);
+    anpipnct->AddCorrFctn(csqqinvpipnct);
+
+    AliFemtoCorrFctnDirectYlm *cylmpipnct = new AliFemtoCorrFctnDirectYlm("cylmpipnct",3,60,0.0,0.3);
+    anpipnct->AddCorrFctn(cylmpipnct);
+    
+    AliFemtoQinvCorrFctn *cqinvpipnct = new AliFemtoQinvCorrFctn("qinvcfpipnct", 100,0.0,1.0);
+    anpipnct->AddCorrFctn(cqinvpipnct);
+
+    AliFemtoChi2CorrFctn *cchiqinvpipnct= new AliFemtoChi2CorrFctn("chicfpipnct",40,0.0,0.4);
+    anpipnct->AddCorrFctn(cchiqinvpipnct);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspipnct = new AliFemtoCorrFctnTPCNcls("cqtpcnclspipnct",40,0.0,0.4);
+    anpipnct->AddCorrFctn(cqtpcnclspipnct);
+
+    Manager->AddAnalysis(anpipnct);    
+  }
+  // *** End pion-pion (positive) analysis
+
+  // *** Begin pion-pion (negative) analysis ***
+  if (runNegativeQA) {
+    AliFemtoVertexMultAnalysis *anpimnct = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpimnct->SetNumEventsToMix(10);
+    anpimnct->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpimnct = new AliFemtoBasicEventCut();
+    mecpimnct->SetEventMult(2,100000);
+    mecpimnct->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpimnct = new AliFemtoESDTrackCut();
+    dtcpimnct->SetPidProbPion(0.0,1.001);
+    dtcpimnct->SetPidProbMuon(0.0,1.0);
+    dtcpimnct->SetPidProbKaon(0.0,1.0);
+    dtcpimnct->SetPidProbProton(0.0,1.0);
+    dtcpimnct->SetCharge(-1.0);
+    dtcpimnct->SetPt(0.1,1.0);
+    dtcpimnct->SetMass(PionMass);
+    // Track quality cuts
+    dtcpimnct->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    //  dtcpimnct->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpimnct->SetminTPCncls(50);
+    dtcpimnct->SetRemoveKinks(kTRUE);
+    dtcpimnct->SetLabel(kFALSE);
+    dtcpimnct->SetMaxITSChiNdof(6.0);
+    dtcpimnct->SetMaxTPCChiNdof(6.0);
+    dtcpimnct->SetMaxImpactXY(3.0);
+    dtcpimnct->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpimnct = new AliFemtoCutMonitorParticleYPt("cutPasspimnct", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpimnct = new AliFemtoCutMonitorParticleYPt("cutFailpimnct", 0.13957);
+    dtcpimnct->AddCutMonitor(cutPassYPtpimnct, cutFailYPtpimnct);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimnct = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpimnct->SetShareQualityMax(1.0);
+    sqpcpimnct->SetShareFractionMax(1.0);
+    sqpcpimnct->SetRemoveSameLabel(kFALSE);
+    sqpcpimnct->SetTPCEntranceSepMinimum(0.0);
+
+    anpimnct->SetEventCut(mecpimnct);
+    anpimnct->SetFirstParticleCut(dtcpimnct);
+    anpimnct->SetSecondParticleCut(dtcpimnct);
+    anpimnct->SetPairCut(sqpcpimnct);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpimnct= new AliFemtoShareQualityCorrFctn("sqqinvcfpimnct",40,0.0,0.4);
+    anpimnct->AddCorrFctn(csqqinvpimnct);
+
+    AliFemtoCorrFctnDirectYlm *cylmpimnct = new AliFemtoCorrFctnDirectYlm("cylmpimnct",3,60,0.0,0.3);
+    anpimnct->AddCorrFctn(cylmpimnct);
+    
+    AliFemtoQinvCorrFctn *cqinvpimnct = new AliFemtoQinvCorrFctn("qinvcfpimnct", 100,0.0,1.0);
+    anpimnct->AddCorrFctn(cqinvpimnct);
+
+    AliFemtoChi2CorrFctn *cchiqinvpimnct= new AliFemtoChi2CorrFctn("chicfpimnct",40,0.0,0.4);
+    anpimnct->AddCorrFctn(cchiqinvpimnct);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspimnct = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimnct",40,0.0,0.4);
+    anpimnct->AddCorrFctn(cqtpcnclspimnct);
+
+    Manager->AddAnalysis(anpimnct);    
+  }
+  // *** End pion-pion (negative) analysis
+
+  // *** Third QA task - HBT analysis with all pair cuts off, TPC only ***
+  // *** Begin pion-pion (positive) analysis ***
+  if (runPositiveTPCQA) {
+    AliFemtoVertexMultAnalysis *anpiptpc = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpiptpc->SetNumEventsToMix(10);
+    anpiptpc->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpiptpc = new AliFemtoBasicEventCut();
+    mecpiptpc->SetEventMult(2,100000);
+    mecpiptpc->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpiptpc = new AliFemtoESDTrackCut();
+    dtcpiptpc->SetPidProbPion(0.0,1.001);
+    dtcpiptpc->SetPidProbMuon(0.0,1.0);
+    dtcpiptpc->SetPidProbKaon(0.0,1.0);
+    dtcpiptpc->SetPidProbProton(0.0,1.0);
+    dtcpiptpc->SetCharge(1.0);
+    dtcpiptpc->SetPt(0.1,1.0);
+    dtcpiptpc->SetMass(PionMass);
+    // Track quality cuts
+    //dtcpiptpc->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    dtcpiptpc->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpiptpc->SetminTPCncls(50);
+    dtcpiptpc->SetRemoveKinks(kTRUE);
+    dtcpiptpc->SetLabel(kFALSE);
+    dtcpiptpc->SetMaxITSChiNdof(6.0);
+    dtcpiptpc->SetMaxTPCChiNdof(6.0);
+    dtcpiptpc->SetMaxImpactXY(3.0);
+    dtcpiptpc->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpiptpc = new AliFemtoCutMonitorParticleYPt("cutPasspiptpc", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpiptpc = new AliFemtoCutMonitorParticleYPt("cutFailpiptpc", 0.13957);
+    dtcpiptpc->AddCutMonitor(cutPassYPtpiptpc, cutFailYPtpiptpc);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpiptpc = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpiptpc->SetShareQualityMax(1.0);
+    sqpcpiptpc->SetShareFractionMax(1.0);
+    sqpcpiptpc->SetRemoveSameLabel(kFALSE);
+    sqpcpiptpc->SetTPCEntranceSepMinimum(0.0);
+
+    anpiptpc->SetEventCut(mecpiptpc);
+    anpiptpc->SetFirstParticleCut(dtcpiptpc);
+    anpiptpc->SetSecondParticleCut(dtcpiptpc);
+    anpiptpc->SetPairCut(sqpcpiptpc);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpiptpc= new AliFemtoShareQualityCorrFctn("sqqinvcfpiptpc",40,0.0,0.4);
+    anpiptpc->AddCorrFctn(csqqinvpiptpc);
+
+    AliFemtoCorrFctnDirectYlm *cylmpiptpc = new AliFemtoCorrFctnDirectYlm("cylmpiptpc",3,60,0.0,0.3);
+    anpiptpc->AddCorrFctn(cylmpiptpc);
+    
+    AliFemtoQinvCorrFctn *cqinvpiptpc = new AliFemtoQinvCorrFctn("qinvcfpiptpc", 100,0.0,1.0);
+    anpiptpc->AddCorrFctn(cqinvpiptpc);
+
+    AliFemtoChi2CorrFctn *cchiqinvpiptpc= new AliFemtoChi2CorrFctn("chicfpiptpc",40,0.0,0.4);
+    anpiptpc->AddCorrFctn(cchiqinvpiptpc);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspiptpc = new AliFemtoCorrFctnTPCNcls("cqtpcnclspiptpc",40,0.0,0.4);
+    anpiptpc->AddCorrFctn(cqtpcnclspiptpc);
+
+    Manager->AddAnalysis(anpiptpc);    
+  }
+  // *** End pion-pion (positive) analysis
+
+  // *** Begin pion-pion (negative) analysis ***
+  if (runNegativeTPCQA) {
+    AliFemtoVertexMultAnalysis *anpimtpc = new AliFemtoVertexMultAnalysis(3, -15.6, 15.6, 5, 2, 200);
+    anpimtpc->SetNumEventsToMix(10);
+    anpimtpc->SetMinSizePartCollection(2);
+
+    AliFemtoBasicEventCut* mecpimtpc = new AliFemtoBasicEventCut();
+    mecpimtpc->SetEventMult(2,100000);
+    mecpimtpc->SetVertZPos(-1000,1000);
+       
+    AliFemtoESDTrackCut* dtcpimtpc = new AliFemtoESDTrackCut();
+    dtcpimtpc->SetPidProbPion(0.0,1.001);
+    dtcpimtpc->SetPidProbMuon(0.0,1.0);
+    dtcpimtpc->SetPidProbKaon(0.0,1.0);
+    dtcpimtpc->SetPidProbProton(0.0,1.0);
+    dtcpimtpc->SetCharge(-1.0);
+    dtcpimtpc->SetPt(0.1,1.0);
+    dtcpimtpc->SetMass(PionMass);
+    // Track quality cuts
+    //dtcpimtpc->SetStatus(AliESDtrack::kTPCrefit|AliESDtrack::kITSrefit);
+    dtcpimtpc->SetStatus(AliESDtrack::kTPCrefit);
+    dtcpimtpc->SetminTPCncls(50);
+    dtcpimtpc->SetRemoveKinks(kTRUE);
+    dtcpimtpc->SetLabel(kFALSE);
+    dtcpimtpc->SetMaxITSChiNdof(6.0);
+    dtcpimtpc->SetMaxTPCChiNdof(6.0);
+    dtcpimtpc->SetMaxImpactXY(3.0);
+    dtcpimtpc->SetMaxImpactZ(3.0);
+
+    AliFemtoCutMonitorParticleYPt *cutPassYPtpimtpc = new AliFemtoCutMonitorParticleYPt("cutPasspimtpc", 0.13957);
+    AliFemtoCutMonitorParticleYPt *cutFailYPtpimtpc = new AliFemtoCutMonitorParticleYPt("cutFailpimtpc", 0.13957);
+    dtcpimtpc->AddCutMonitor(cutPassYPtpimtpc, cutFailYPtpimtpc);
+
+    AliFemtoShareQualityTPCEntranceSepPairCut *sqpcpimtpc = new AliFemtoShareQualityTPCEntranceSepPairCut();
+    sqpcpimtpc->SetShareQualityMax(1.0);
+    sqpcpimtpc->SetShareFractionMax(1.0);
+    sqpcpimtpc->SetRemoveSameLabel(kFALSE);
+    sqpcpimtpc->SetTPCEntranceSepMinimum(0.0);
+
+    anpimtpc->SetEventCut(mecpimtpc);
+    anpimtpc->SetFirstParticleCut(dtcpimtpc);
+    anpimtpc->SetSecondParticleCut(dtcpimtpc);
+    anpimtpc->SetPairCut(sqpcpimtpc);
+    
+    AliFemtoShareQualityCorrFctn *csqqinvpimtpc= new AliFemtoShareQualityCorrFctn("sqqinvcfpimtpc",40,0.0,0.4);
+    anpimtpc->AddCorrFctn(csqqinvpimtpc);
+
+    AliFemtoCorrFctnDirectYlm *cylmpimtpc = new AliFemtoCorrFctnDirectYlm("cylmpimtpc",3,60,0.0,0.3);
+    anpimtpc->AddCorrFctn(cylmpimtpc);
+    
+    AliFemtoQinvCorrFctn *cqinvpimtpc = new AliFemtoQinvCorrFctn("qinvcfpimtpc", 100,0.0,1.0);
+    anpimtpc->AddCorrFctn(cqinvpimtpc);
+
+    AliFemtoChi2CorrFctn *cchiqinvpimtpc= new AliFemtoChi2CorrFctn("chicfpimtpc",40,0.0,0.4);
+    anpimtpc->AddCorrFctn(cchiqinvpimtpc);
+
+    AliFemtoCorrFctnTPCNcls *cqtpcnclspimtpc = new AliFemtoCorrFctnTPCNcls("cqtpcnclspimtpc",40,0.0,0.4);
+    anpimtpc->AddCorrFctn(cqtpcnclspimtpc);
+
+    Manager->AddAnalysis(anpimtpc);    
+  }
+  // *** End pion-pion (negative) analysis
+
+
+
+  return Manager;
+}                         
+                      
diff --git a/PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/FemtoQAPlots.C b/PWG2/FEMTOSCOPY/macros/Train/TwoTrackQA/FemtoQAPlots.C
new file mode 100644 (file)
index 0000000..7ac92f9
--- /dev/null
@@ -0,0 +1,423 @@
+void preparepad()
+{
+  gPad->SetFillColor(0);
+  gPad->SetFillStyle(4000);
+
+  gPad->SetTopMargin(0.02);
+  //  gPad->SetRightMargin(0.02);
+}
+
+void preparehist(TH1D *hist)
+{
+  hist->SetMarkerSize(0.8);
+  hist->SetMarkerStyle(8);
+  hist->SetMarkerColor(2);
+}
+
+TH2D *plot2ddtpc(TH2D *numdtpc, TH2D *dendtpc, const char *suff)
+{
+  // Create a 2D CF vs minimum separation distance in the TPC
+
+
+
+  int nslices = numdtpc->GetNbinsY();
+  int nrange = numdtpc->GetNbinsX();
+  char hname[200];
+  
+  sprintf(hname,"numdtpc%s",suff);
+  TH2D *nums = new TH2D(hname, ";q_{inv} [GeV/c];nominal TPC entrance sep - lower cut-off [cm]", 
+                       numdtpc->GetNbinsX(), numdtpc->GetXaxis()->GetXmin(), numdtpc->GetXaxis()->GetXmax(), 
+                       numdtpc->GetNbinsY(), numdtpc->GetYaxis()->GetXmin(), numdtpc->GetYaxis()->GetXmax());
+  sprintf(hname,"dendtpc%s",suff);
+  TH2D *dens = new TH2D(hname, ";q_{inv} [GeV/c];nominal TPC entrance sep - lower cut-off [cm]", 
+                       numdtpc->GetNbinsX(), numdtpc->GetXaxis()->GetXmin(), numdtpc->GetXaxis()->GetXmax(), 
+                       numdtpc->GetNbinsY(), numdtpc->GetYaxis()->GetXmin(), numdtpc->GetYaxis()->GetXmax());
+  sprintf(hname,"ratdtpc%s",suff);
+  TH2D *rats = new TH2D(hname,  ";q_{inv} [GeV/c];nominal TPC entrance sep - lower cut-off [cm]", 
+                       numdtpc->GetNbinsX(), numdtpc->GetXaxis()->GetXmin(), numdtpc->GetXaxis()->GetXmax(), 
+                       numdtpc->GetNbinsY(), numdtpc->GetYaxis()->GetXmin(), numdtpc->GetYaxis()->GetXmax());
+  
+  sprintf(hname,"chi2toone%s",suff);
+  TH1D *chi2toone = new TH1D(hname,";nominal TPC entrance sep - lower cut-off [cm];\\chi^2",
+                            numdtpc->GetNbinsY(), numdtpc->GetYaxis()->GetXmin(), numdtpc->GetYaxis()->GetXmax());
+
+  char bufname[100];
+  for (int iter=0; iter<nslices; iter++) {
+    sprintf(bufname, "numbuf%i", iter);
+    TH1D *numbuf = numdtpc->ProjectionX(bufname, (iter+1),  nslices, "e");
+    sprintf(bufname, "denbuf%i", iter);
+    TH1D *denbuf = dendtpc->ProjectionX(bufname, (iter+1),  nslices, "e");
+    
+    TH1D *ratbuf = new TH1D(*numbuf);
+    ratbuf->Divide(denbuf);
+    Double_t scale = ratbuf->Integral(nrange-9,nrange)/10.0;
+
+    if (scale > 0.0)
+      ratbuf->Scale(1.0/scale);
+    
+    double chi2 = 0.0;
+    for (int ibin = 1; ibin<=ratbuf->GetNbinsX(); ibin++) {
+      rats->SetBinContent(ibin, iter+1, ratbuf->GetBinContent(ibin));
+      if (ratbuf->GetBinError(ibin) > 0.0) {
+       chi2 += (TMath::Power(ratbuf->GetBinContent(ibin)-1.0,2)/
+                TMath::Power(ratbuf->GetBinError(ibin),2));
+      }
+    }
+    chi2toone->SetBinContent(iter+1, chi2);
+  }
+
+  rats->GetXaxis()->SetTitleOffset(0.9);
+  rats->GetYaxis()->SetTitleOffset(0.95);
+
+  return rats;
+}
+
+TH2D * plot2dshare(TH2D *numshare, TH2D *denshare, const char *suff)
+{
+  // Create a 2D CF vs minimum separation distance in the TPC
+  int nrange = numshare->GetNbinsX();
+  char hname[200];
+    
+  sprintf(hname,"numshare%s",suff);
+  TH2D *nums = new TH2D(hname, ";q_{inv} [GeV/c];fraction of shared hits - upper cut-off", 
+                       numshare->GetNbinsX(), numshare->GetXaxis()->GetXmin(), numshare->GetXaxis()->GetXmax(),
+                       numshare->GetNbinsY(), numshare->GetYaxis()->GetXmin(), numshare->GetYaxis()->GetXmax());
+
+  sprintf(hname,"denshare%s",suff);
+  TH2D *dens = new TH2D(hname, ";q_{inv} [GeV/c];fraction of shared hits - upper cut-off", 
+                       numshare->GetNbinsX(), numshare->GetXaxis()->GetXmin(), numshare->GetXaxis()->GetXmax(),
+                       numshare->GetNbinsY(), numshare->GetYaxis()->GetXmin(), numshare->GetYaxis()->GetXmax());
+
+  sprintf(hname,"ratshare%s",suff);
+  TH2D *rats = new TH2D(hname, ";q_{inv} [GeV/c];fraction of shared hits - upper cut-off", 
+                       numshare->GetNbinsX(), numshare->GetXaxis()->GetXmin(), numshare->GetXaxis()->GetXmax(),
+                       numshare->GetNbinsY(), numshare->GetYaxis()->GetXmin(), numshare->GetYaxis()->GetXmax());
+  
+  char bufname[100];
+  for (int iter=0; iter<numshare->GetNbinsY(); iter++) {
+    sprintf(bufname, "numbuf%i", iter);
+    TH1D *numbuf = numshare->ProjectionX(bufname, 1, iter+1, "e");
+    sprintf(bufname, "denbuf%i", iter);
+    TH1D *denbuf = denshare->ProjectionX(bufname, 1, iter+1, "e");
+    
+    ratbuf = new TH1D(*numbuf);
+    ratbuf->Divide(denbuf);
+    Double_t scale = ratbuf->Integral(nrange-24,nrange)/25.0;
+    ratbuf->Scale(1.0/scale);
+
+    for (int ibin = 1; ibin<=ratbuf->GetNbinsX(); ibin++) {
+      rats->SetBinContent(ibin, iter+1, ratbuf->GetBinContent(ibin));
+    }
+  }
+  
+  rats->GetXaxis()->SetTitleOffset(0.95);
+  rats->GetYaxis()->SetTitleOffset(0.9);
+
+  return rats;
+}
+
+TH2D *plot2dquality(TH2D *numqual, TH2D *denqual, const char *suff)
+{
+  // Create a 2D CF vs pair quality 
+  int nrange = numqual->GetNbinsX();
+  char hname[200];
+    
+  sprintf(hname,"numqual%s",suff);
+  TH2D *nums = new TH2D(hname, ";q_{inv} [GeV/c];pair quality - upper cut-off", 
+                       numqual->GetNbinsX(), numqual->GetXaxis()->GetXmin(), numqual->GetXaxis()->GetXmax(),
+                       numqual->GetNbinsY(), numqual->GetYaxis()->GetXmin(), numqual->GetYaxis()->GetXmax());
+  sprintf(hname,"denqual%s",suff);
+  TH2D *dens = new TH2D(hname, ";q_{inv} [GeV/c];pair quality - upper cut-off", 
+                       numqual->GetNbinsX(), numqual->GetXaxis()->GetXmin(), numqual->GetXaxis()->GetXmax(),
+                       numqual->GetNbinsY(), numqual->GetYaxis()->GetXmin(), numqual->GetYaxis()->GetXmax());
+  sprintf(hname,"ratqual%s",suff);
+  TH2D *rats = new TH2D(hname, ";q_{inv} [GeV/c];pair quality - upper cut-off", 
+                       numqual->GetNbinsX(), numqual->GetXaxis()->GetXmin(), numqual->GetXaxis()->GetXmax(),
+                       numqual->GetNbinsY(), numqual->GetYaxis()->GetXmin(), numqual->GetYaxis()->GetXmax());
+  
+  char bufname[100];
+  for (int iter=0; iter<numqual->GetNbinsY(); iter++) {
+    sprintf(bufname, "numbuf%i", iter);
+    TH1D *numbuf = numqual->ProjectionX(bufname, 1, iter+1, "e");
+    sprintf(bufname, "denbuf%i", iter);
+    TH1D *denbuf = denqual->ProjectionX(bufname, 1, iter+1, "e");
+    
+    ratbuf = new TH1D(*numbuf);
+    ratbuf->Divide(denbuf);
+    Double_t scale = ratbuf->Integral(nrange-24,nrange)/25.0;
+    if (scale > 0.0)
+      ratbuf->Scale(1.0/scale);
+
+    for (int ibin = 1; ibin<=ratbuf->GetNbinsX(); ibin++) {
+      rats->SetBinContent(ibin, iter+1, ratbuf->GetBinContent(ibin));
+    }
+  }
+  
+  rats->GetXaxis()->SetTitleOffset(0.95);
+  rats->GetYaxis()->SetTitleOffset(0.9);
+
+  return rats;
+}
+
+TH1D *getcf(TH1D *numq, TH1D *denq, const char *suff)
+{
+  char hname[200];
+
+  TH1D *rats = new TH1D(*numq);
+  rats->Reset("ICE");
+  rats->Divide(numq, denq, 1.0*denq->Integral()/numq->Integral(), 1.0);
+
+  sprintf(hname, "ratqinv%s", suff);
+  rats->SetName(hname);
+  rats->SetTitle(";q_{inv} [GeV/c];C(q_{inv})");
+  preparehist(rats);
+
+  return rats;
+}
+
+void FemtoQAPlots(const char *infname)
+{
+  gStyle->SetStatX(0.8);
+  gStyle->SetStatW(0.25);
+  gStyle->SetOptStat(11);
+
+  TFile *ifile = new TFile(infname);
+  
+  gStyle->SetPalette(1,0);
+
+  // Make plot for TPC entrance separation - positive
+
+  TCanvas *candtpcp = new TCanvas("candtpcp","candtpcp",1200,300);
+  preparepad();
+  candtpcp->Divide(3,1,0.0001,0.0001);
+  
+  candtpcp->cd(1);  preparepad();
+  TH2D *ratdtpcpipstd = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpipstd"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpipstd"), "pipstd");
+  ratdtpcpipstd->Draw("COLZ");
+
+  candtpcp->cd(2);  preparepad();
+  TH2D *ratdtpcpipnct = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpipnct"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpipnct"), "pipnct");
+  ratdtpcpipnct->Draw("COLZ");
+
+  candtpcp->cd(3);  preparepad();
+  TH2D *ratdtpcpiptpc = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpiptpc"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpiptpc"), "piptpc");
+  ratdtpcpiptpc->Draw("COLZ");
+
+  // Make plot for TPC entrance separation - negative
+
+  TCanvas *candtpcm = new TCanvas("candtpcm","candtpcm",1200,300);
+  preparepad();
+  candtpcm->Divide(3,1,0.0001,0.0001);
+  
+  candtpcm->cd(1);  preparepad();
+  TH2D *ratdtpcpimstd = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpimstd"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpimstd"), "pimstd");
+  ratdtpcpimstd->Draw("COLZ");
+
+  candtpcm->cd(2);  preparepad();
+  TH2D *ratdtpcpimnct = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpimnct"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpimnct"), "pimnct");
+  ratdtpcpimnct->Draw("COLZ");
+
+  candtpcm->cd(3);  preparepad();
+  TH2D *ratdtpcpimtpc = plot2ddtpc((TH2D *) gDirectory->Get("NumTPCSepsqqinvcfpimtpc"), (TH2D *) gDirectory->Get("DenTPCSepsqqinvcfpimtpc"), "pimtpc");
+  ratdtpcpimtpc->Draw("COLZ");
+
+  // Make plot for pair fraction of shared hits - positive
+
+  TCanvas *cansharep = new TCanvas("cansharep","cansharep",1200,300);
+  preparepad();
+  cansharep->Divide(3,1,0.0001,0.0001);
+  
+  cansharep->cd(1);  preparepad();
+  TH2D *ratsharepipstd = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpipstd"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpipstd"), "pipstd");
+  ratsharepipstd->Draw("COLZ");
+
+  cansharep->cd(2);  preparepad();
+  TH2D *ratsharepipnct = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpipnct"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpipnct"), "pipnct");
+  ratsharepipnct->Draw("COLZ");
+
+  cansharep->cd(3);  preparepad();
+  TH2D *ratsharepiptpc = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpiptpc"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpiptpc"), "piptpc");
+  ratsharepiptpc->Draw("COLZ");
+
+  // Make plot for pair fraction of shared hits - negative
+
+  TCanvas *cansharem = new TCanvas("cansharem","cansharem",1200,300);
+  preparepad();
+  cansharem->Divide(3,1,0.0001,0.0001);
+  
+  cansharem->cd(1);  preparepad();
+  TH2D *ratsharepimstd = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpimstd"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpimstd"), "pimstd");
+  ratsharepimstd->Draw("COLZ");
+
+  cansharem->cd(2);  preparepad();
+  TH2D *ratsharepimnct = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpimnct"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpimnct"), "pimnct");
+  ratsharepimnct->Draw("COLZ");
+
+  cansharem->cd(3);  preparepad();
+  TH2D *ratsharepimtpc = plot2dshare((TH2D *) gDirectory->Get("NumSharesqqinvcfpimtpc"), (TH2D *) gDirectory->Get("DenSharesqqinvcfpimtpc"), "pimtpc");
+  ratsharepimtpc->Draw("COLZ");
+
+  // Make signal pair fraction of shared hits - positive
+
+  TCanvas *cansharesignalp = new TCanvas("cansharesignalp","cansharesignalp",1200,300);
+  preparepad();
+  cansharesignalp->Divide(3,1,0.0001,0.0001);
+  
+  cansharesignalp->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpipstd"))->Draw("COLZ");
+
+  cansharesignalp->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpipnct"))->Draw("COLZ");
+
+  cansharesignalp->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpiptpc"))->Draw("COLZ");
+
+  // Make signal pair fraction of shared hits - negative
+
+  TCanvas *cansharesignalm = new TCanvas("cansharesignalm","cansharesignalm",1200,300);
+  preparepad();
+  cansharesignalm->Divide(3,1,0.0001,0.0001);
+  
+  cansharesignalm->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpimstd"))->Draw("COLZ");
+
+  cansharesignalm->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpimnct"))->Draw("COLZ");
+
+  cansharesignalm->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("NumSharesqqinvcfpimtpc"))->Draw("COLZ");
+
+  // Make background pair fraction of shared hits - positive
+
+  TCanvas *cansharebackp = new TCanvas("cansharebackp","cansharebackp",1200,300);
+  preparepad();
+  cansharebackp->Divide(3,1,0.0001,0.0001);
+  
+  cansharebackp->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpipstd"))->Draw("COLZ");
+
+  cansharebackp->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpipnct"))->Draw("COLZ");
+
+  cansharebackp->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpiptpc"))->Draw("COLZ");
+
+  // Make background pair fraction of shared hits - negative
+
+  TCanvas *cansharebackm = new TCanvas("cansharebackm","cansharebackm",1200,300);
+  preparepad();
+  cansharebackm->Divide(3,1,0.0001,0.0001);
+  
+  cansharebackm->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpimstd"))->Draw("COLZ");
+
+  cansharebackm->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpimnct"))->Draw("COLZ");
+
+  cansharebackm->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("DenSharesqqinvcfpimtpc"))->Draw("COLZ");
+
+  // Make plot for pair quality - positive
+
+  TCanvas *canqualp = new TCanvas("canqualp","canqualp",1200,300);
+  preparepad();
+  canqualp->Divide(3,1,0.0001,0.0001);
+  
+  canqualp->cd(1);  preparepad();
+  TH2D *ratqualpipstd = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpipstd"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpipstd"), "pipstd");
+  ratqualpipstd->Draw("COLZ");
+
+  canqualp->cd(2);  preparepad();
+  TH2D *ratqualpipnct = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpipnct"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpipnct"), "pipnct");
+  ratqualpipnct->Draw("COLZ");
+
+  canqualp->cd(3);  preparepad();
+  TH2D *ratqualpiptpc = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpiptpc"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpiptpc"), "piptpc");
+  ratqualpiptpc->Draw("COLZ");
+
+  // Make plot for pair quality - negative
+
+  TCanvas *canqualm = new TCanvas("canqualm","canqualm",1200,300);
+  preparepad();
+  canqualm->Divide(3,1,0.0001,0.0001);
+  
+  canqualm->cd(1);  preparepad();
+  TH2D *ratqualpimstd = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimstd"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpimstd"), "pimstd");
+  ratqualpimstd->Draw("COLZ");
+
+  canqualm->cd(2);  preparepad();
+  TH2D *ratqualpimnct = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimnct"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpimnct"), "pimnct");
+  ratqualpimnct->Draw("COLZ");
+
+  canqualm->cd(3);  preparepad();
+  TH2D *ratqualpimtpc = plot2dquality((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimtpc"), (TH2D *) gDirectory->Get("DenQualitysqqinvcfpimtpc"), "pimtpc");
+  ratqualpimtpc->Draw("COLZ");
+
+  // Make signal pair quality - positive
+
+  TCanvas *canqualitysignalp = new TCanvas("canqualitysignalp","canqualitysignalp",1200,300);
+  preparepad();
+  canqualitysignalp->Divide(3,1,0.0001,0.0001);
+  
+  canqualitysignalp->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpipstd"))->Draw("COLZ");
+
+  canqualitysignalp->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpipnct"))->Draw("COLZ");
+
+  canqualitysignalp->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpiptpc"))->Draw("COLZ");
+
+  // Make signal pair quality - negative
+
+  TCanvas *canqualitysignalm = new TCanvas("canqualitysignalm","canqualitysignalm",1200,300);
+  preparepad();
+  canqualitysignalm->Divide(3,1,0.0001,0.0001);
+  
+  canqualitysignalm->cd(1);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimstd"))->Draw("COLZ");
+
+  canqualitysignalm->cd(2);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimnct"))->Draw("COLZ");
+
+  canqualitysignalm->cd(3);  preparepad();
+  ((TH2D *) gDirectory->Get("NumQualitysqqinvcfpimtpc"))->Draw("COLZ");
+
+  // Make plot for pair quality - positive
+
+  TCanvas *canqinvp = new TCanvas("canqinvp","canqinvp",1200,300);
+  preparepad();
+  canqinvp->Divide(3,1,0.0001,0.0001);
+  
+  canqinvp->cd(1); preparepad();
+  TH1D *ratqinvpipstd = getcf((TH1D *) gDirectory->Get("Numqinvcfpipstd"), (TH1D *) gDirectory->Get("Denqinvcfpipstd"), "pipstd");
+  ratqinvpipstd->Draw();
+
+  canqinvp->cd(2); preparepad();
+  TH1D *ratqinvpipnct = getcf((TH1D *) gDirectory->Get("Numqinvcfpipnct"), (TH1D *) gDirectory->Get("Denqinvcfpipnct"), "pipnct");
+  ratqinvpipnct->Draw("");
+
+  canqinvp->cd(3); preparepad();
+  TH1D *ratqinvpiptpc = getcf((TH1D *) gDirectory->Get("Numqinvcfpiptpc"), (TH1D *) gDirectory->Get("Denqinvcfpiptpc"), "piptpc");
+  ratqinvpiptpc->Draw("");
+
+  // Make plot for pair quality - negative
+
+  TCanvas *canqinvm = new TCanvas("canqinvm","canqinvm",1200,300);
+  preparepad();
+  canqinvm->Divide(3,1,0.0001,0.0001);
+  
+  canqinvm->cd(1); preparepad();
+  TH1D *ratqinvpimstd = getcf((TH1D *) gDirectory->Get("Numqinvcfpimstd"), (TH1D *) gDirectory->Get("Denqinvcfpimstd"), "pimstd");
+  ratqinvpimstd->Draw("");
+
+  canqinvm->cd(2); preparepad();
+  TH1D *ratqinvpimnct = getcf((TH1D *) gDirectory->Get("Numqinvcfpimnct"), (TH1D *) gDirectory->Get("Denqinvcfpimnct"), "pimnct");
+  ratqinvpimnct->Draw("");
+
+  canqinvm->cd(3); preparepad();
+  TH1D *ratqinvpimtpc = getcf((TH1D *) gDirectory->Get("Numqinvcfpimtpc"), (TH1D *) gDirectory->Get("Denqinvcfpimtpc"), "pimtpc");
+  ratqinvpimtpc->Draw("");
+
+  
+}