new attempt with eventpoolmanager. fixes strange branch problem.
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Sep 2010 10:00:37 +0000 (10:00 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 24 Sep 2010 10:00:37 +0000 (10:00 +0000)
PWG4/TwoPartCorr/anaCorr.C

index a404493..0906d7a 100644 (file)
 #include <TH2F.h>
 #include <TTimeStamp.h>
 #include "TreeClasses.h"
-#include "EventPool.h"
+#include "EventPoolManager.h"
 #endif
 
-void anaCorr(Int_t nEvents=-1,
-               Int_t nmin=100,
-               Int_t nmax=500,
-               Double_t vzmin=-10,
-               Double_t vzmax=+10,
-               Double_t ptmin1=0.150,
-               Double_t ptmax1=10,
-               Double_t ptmin2=0.150,
-               Double_t ptmax2=10,
-               Double_t etamin1=-1.4,
-               Double_t etamax1=+1.4,
-               Double_t etamin2=-1.4,
-               Double_t etamax2=+1.4,
-               Double_t psize=100,
-               Double_t nmix=10,
-               const char *inFileNames = "res/mergedruns/merged_run*.root");
+void anaCorr(const char *inFileNames,
+             Int_t    nEvents=-1,
+             Int_t    nmin=10,
+             Int_t    nmax=100,
+             Double_t vzmin=-8,
+             Double_t vzmax=+8,
+             Int_t    vzbins=8,
+             Int_t    vcmin=10,
+             Int_t    psize=100,
+             Double_t nmix=10,
+             Double_t ptmin=0.4, 
+             Double_t ptmax=10, 
+             Double_t etamin=-1.4, 
+             Double_t etamax=+1.4);
 
 inline Double_t DeltaPhi(const MyPart &t1, const MyPart &t2, 
                          Double_t rangeMin = -TMath::Pi()/2, 
@@ -37,27 +35,41 @@ inline Double_t DeltaEta(const MyPart &t1, const MyPart &t2);
 inline Bool_t   InBounds(Double_t val, Double_t min, Double_t max);
 inline Double_t InvMass(const MyPart &p1, const MyPart &p2);
 
+class Noti: public TObject
+{
+public:
+  Noti() : fc(0) {;}
+  virtual ~Noti() {;}
+  Bool_t Notify()         { fc=1; return 1; }
+  Bool_t Notified() const { return fc;      }
+  void   Reset()          { fc=0;           }
+protected:
+  Bool_t fc; //=1 when file changed
+  ClassDef (Noti,0)
+};
+
 //-----------------------------------------------------------------------------------------------------
 
-void anaCorr(Int_t nEvents,
-             Int_t nmin,
-             Int_t nmax,
+void anaCorr(const char *inFileNames,
+             Int_t    nEvents,
+             Int_t    nmin,
+             Int_t    nmax,
              Double_t vzmin,
              Double_t vzmax,
-             Double_t ptmin1,
-             Double_t ptmax1,
-             Double_t ptmin2,
-             Double_t ptmax2,
-             Double_t etamin1,
-             Double_t etamax1,
-             Double_t etamin2,
-             Double_t etamax2,
-             Double_t psize,
+             Int_t    vzbins,
+             Int_t    vcmin,
+             Int_t    psize,
              Double_t nmix,
-             const char *inFileNames)
+             Double_t ptmin,
+             Double_t ptmax,
+             Double_t etamin,
+             Double_t etamax)
 {
+  Noti *nme = new Noti;
   TChain *tree = new TChain("MyTree");
   tree->Add(inFileNames);
+  tree->SetNotify(nme);
+
   Int_t nents = tree->GetEntries();
   cout << "Found " << nents << " entries!" << endl;
   if (nents<=0)
@@ -69,14 +81,29 @@ void anaCorr(Int_t nEvents,
 
   if (TClass::GetClass("MyPart"))
     TClass::GetClass("MyPart")->IgnoreTObjectStreamer();
-  MyHeader   *header = 0;
-  TClonesArray *tracks = 0;
-  TBranch *hBranch = tree->GetBranch("header");
-  hBranch->SetAddress(&header);
-  TBranch* tBranch = tree->GetBranch("parts");
-  tBranch->SetAddress(&tracks);
-
-  EventPool *pool = new EventPool(psize);
+  if (TClass::GetClass("MyTracklet"))
+    TClass::GetClass("MyTracklet")->IgnoreTObjectStreamer();
+
+  MyHeader     *header  = 0;
+  TClonesArray *tracks  = 0;
+  TBranch      *hBranch = 0;
+  TBranch      *tBranch = 0;
+
+  Double_t *multbinsa = new Double_t[2]; multbinsa[0] = nmin;  multbinsa[1] = nmax; 
+  Double_t *vzbinsa = new Double_t[vzbins+1];
+  Double_t vzbinwidth = Double_t(vzmax-vzmin)/vzbins;
+  for(Int_t i=0;i<vzbins+1;++i) {
+    vzbinsa[i] = vzmin+i*vzbinwidth;
+  }
+  TH1F* hMultBins = new TH1F("hMultBins", "Event multiplicity binning", 
+                            1, multbinsa);
+  TH1F* hZvtxBins = new TH1F("hZvtxBins", "Event Z-vertex binning", 
+                            vzbins, vzbinsa);
+
+
+  EventPoolManager *pm = new EventPoolManager;
+  pm->InitEventPools(psize, 1, multbinsa, vzbins, vzbinsa);
+
   TH2F *hSig = new TH2F("hSig",";#Delta#eta;#Delta#phi",60,-3,3,64,-TMath::Pi()/2,3*TMath::Pi()/2);
   TH2F *hBkg = new TH2F("hBgk",";#Delta#eta;#Delta#phi",60,-3,3,64,-TMath::Pi()/2,3*TMath::Pi()/2);
 
@@ -84,38 +111,58 @@ void anaCorr(Int_t nEvents,
   Double_t demax = 0.01;
 
   for (Int_t i=0;i<nents;++i) {
-    hBranch->GetEntry(i);
+    Int_t li = tree->LoadTree(i);
+    if (nme->Notified()) {
+      hBranch = tree->GetBranch("header");
+      hBranch->SetAddress(&header);
+      tBranch = tree->GetBranch("parts");
+      tBranch->SetAddress(&tracks);
+      nme->Reset();
+    }
+
+    hBranch->GetEntry(li);
     if (header->fIsPileupSPD)
       continue;
     if ((header->fVz<vzmin) ||
         (header->fVz>vzmax) )
       continue;
-//    if ((header->fVc<50) ||
-//        (header->fVc>nmax) )
-//      continue;
-    if ((header->fNSelTracks<nmin) ||
-        (header->fNSelTracks>nmax))
+    if ((header->fVc<vcmin) ||
+        (header->fVc>nmax) )
       continue;
-    tBranch->GetEntry(i);
-
-    if (pool->Depth()==psize) {
-      cout << "Working on event " << i << endl;
-
-      Int_t ntracks = tracks->GetEntries();
-      for (Int_t t1=0;t1<ntracks;++t1) {
-        MyPart *part1 = (MyPart*)tracks->At(t1);
-        if (!InBounds(part1->fPt,ptmin1,ptmax1))
-          continue;
-        if (!InBounds(part1->fEta,etamin1,etamax1))
-          continue;
-        for (Int_t t2=t1+1;t2<ntracks;++t2) {
-          MyPart *part2 = (MyPart*)tracks->At(t2);
-          if (!InBounds(part2->fPt,ptmin2,ptmax2))
-            continue;
-          if (!InBounds(part2->fEta,etamin2,etamax2))
-            continue;
-          //if (InvMass(*part1,*part2)<0.001)
-          //  continue;
+
+    tBranch->GetEntry(li);
+    Int_t ntracks = tracks->GetEntries();
+    TClonesArray *seltracks = new TClonesArray("MyPart",ntracks);
+    Int_t nseltracks = 0;
+    for (Int_t t1=0;t1<ntracks;++t1) {
+      MyPart *part1 = (MyPart*)tracks->At(t1);
+      if (!InBounds(part1->fPt,ptmin,ptmax))
+        continue;
+      if (!InBounds(part1->fEta,etamin,etamax))
+        continue;
+      if (TMath::Abs(part1->fD)>3)
+        continue;
+      if (TMath::Abs(part1->fZ)>3)
+        continue;
+      new((*seltracks)[nseltracks]) MyPart(*part1);
+      ++nseltracks;
+    }
+
+    Int_t mBin = hMultBins->FindBin(nseltracks) - 1;
+    Int_t zBin = hZvtxBins->FindBin(header->fVz) - 1;
+    GenericEventPool *pool = pm->GetEventPool(mBin, zBin);
+    if (!pool) {
+      continue;
+    }
+
+    if (pool->IsReady()) {
+      cout << "Working on event " << i << " from " << header->fRun << endl;
+
+      for (Int_t t1=0;t1<nseltracks;++t1) {
+        MyPart *part1 = (MyPart*)seltracks->At(t1);
+        for (Int_t t2=t1+1;t2<nseltracks;++t2) {
+          MyPart *part2 = (MyPart*)seltracks->At(t2);
+
           Double_t dphi = DeltaPhi(*part1,*part2);
           Double_t deta = DeltaEta(*part1,*part2);
           Double_t dr = dphi*dphi/(dpmax*dpmax) + deta*deta/(demax*demax);
@@ -124,13 +171,7 @@ void anaCorr(Int_t nEvents,
 
           Int_t mtracks = nmix;
           for (Int_t t3=0;t3<mtracks;) {
-            MyPart *part2 = pool->GetRandomTrack();
-            if (!InBounds(part2->fPt,ptmin2,ptmax2))
-              continue;
-            if (!InBounds(part2->fEta,etamin2,etamax2))
-              continue;
-            //if (InvMass(*part1,*part2)<0.001)
-            //  continue;
+            MyPart *part2 = (MyPart*)pool->GetRandomTrack();
             Double_t dphi = DeltaPhi(*part1,*part2);
             Double_t deta = DeltaEta(*part1,*part2);
             Double_t dr = dphi*dphi/(dpmax*dpmax) + deta*deta/(demax*demax);
@@ -140,10 +181,10 @@ void anaCorr(Int_t nEvents,
           }
         }
       }
-    } else {
-      pool->PrintInfo();
-      pool->UpdatePool(i,header,tracks);
     }
+    pool->UpdatePool(i,header,seltracks);
+    if (!pool->IsReady())
+      pool->PrintInfo();
   }
 
   TCanvas *c1 = new TCanvas("cSig");
@@ -183,7 +224,7 @@ Double_t DeltaPhi(const MyPart &t1, const MyPart &t2,
 {
   Double_t dphi = -999;
   Double_t pi = TMath::Pi();
-    Double_t phia = t1.fPhi;  
+  Double_t phia = t1.fPhi;  
   Double_t phib = t2.fPhi;  
   
   if (phia < 0)         phia += 2*pi;