]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Added function to extract good particles
authorvestbo <vestbo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 25 Oct 2001 15:49:09 +0000 (15:49 +0000)
committervestbo <vestbo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 25 Oct 2001 15:49:09 +0000 (15:49 +0000)
HLT/hough/GetGoodParticles.cxx [new file with mode: 0644]
HLT/hough/GetGoodParticles.h [new file with mode: 0644]

diff --git a/HLT/hough/GetGoodParticles.cxx b/HLT/hough/GetGoodParticles.cxx
new file mode 100644 (file)
index 0000000..426757e
--- /dev/null
@@ -0,0 +1,162 @@
+#ifndef __CINT__
+#include "alles.h"
+#include "GetGoodParticles.h"
+#include "AliL3Transform.h"
+#endif
+
+
+void GetGoodParticles(Int_t slice,char *eventfile,char *digitfile)
+{
+
+  Int_t good_number = 70;
+  
+  struct GoodTrack goodtracks[15000];
+  Int_t nt=0;
+  char *fname = "good_tracks";
+  
+  ifstream in(fname);
+  if(in)
+    {
+      printf("Delete you old good_tracks file\n");
+      return;
+    }
+  
+  TFile *evfile = TFile::Open(eventfile);
+  if(!evfile)
+    {
+      printf("No eventfile\n");
+      return;
+    }
+  gAlice = (AliRun*)evfile->Get("gAlice");
+  gAlice->GetEvent(0);
+
+  AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC");
+  AliTPCParam *param = (AliTPCParam*)evfile->Get("75x40_100x60");
+  
+  TPC->SetParam(param);
+  
+  Int_t ver = TPC->IsVersion();
+  if(ver!=2)
+    {
+      printf("Wrong TPC version\n");
+      return;
+    }
+  
+  Int_t zero=TPC->GetParam()->GetZeroSup();
+  
+  Int_t np = gAlice->GetNtrack();
+  Int_t *good = new Int_t[np];
+  for(Int_t ii=0; ii<np; ii++)
+    good[ii] = 0;
+
+  AliL3Transform *transform = new AliL3Transform();
+
+  TFile *digfile = TFile::Open(digitfile);
+  TTree *tree = (TTree*)digfile->Get("TreeD_75x40_100x60_0");
+  AliSimDigits da, *digits=&da;
+  tree->GetBranch("Segment")->SetAddress(&digits);
+  
+  Int_t *count = new Int_t[np]; //np number of particles.
+  Int_t i;
+  for (i=0; i<np; i++) count[i]=0;
+  for(i=0; i<tree->GetEntries(); i++)
+    {
+      if (!tree->GetEvent(i)) continue;
+      Int_t sec,row;
+      param->AdjustSectorRow(digits->GetID(),sec,row);
+      Int_t sl,padrow;
+      transform->Sector2Slice(sl,padrow,sec,row);
+      if(sl < slice) continue;
+      if(sl != slice) break;
+      digits->First();
+      do {
+       Int_t it=digits->CurrentRow(), ip=digits->CurrentColumn();
+       Short_t dig = digits->GetDigit(it,ip);
+       Int_t idx0=digits->GetTrackID(it,ip,0); 
+       Int_t idx1=digits->GetTrackID(it,ip,1);
+       Int_t idx2=digits->GetTrackID(it,ip,2);
+       
+       if (idx0>=0 && dig>=zero) count[idx0]+=1;
+       if (idx1>=0 && dig>=zero) count[idx1]+=1;
+       if (idx2>=0 && dig>=zero) count[idx2]+=1;
+         } while (digits->Next());
+      
+      for (Int_t j=0; j<np; j++) 
+       {
+         if (count[j]>1) //at least two digits at this padrow 
+           good[j]++;
+         
+         count[j]=0;
+       }
+    }
+    
+  delete[] count;
+  
+  
+  //The following code has been taken from offlinemacro->AliTPCComparison.C
+  
+  TTree *TH=gAlice->TreeH();
+  Int_t npart=(Int_t)TH->GetEntries();
+    
+
+
+  while (npart--) {
+    AliTPChit *hit0=0;
+    
+    TPC->ResetHits();
+    TH->GetEvent(npart);
+    AliTPChit * hit = (AliTPChit*) TPC->FirstHit(-1);
+    while (hit){
+      if (hit->fQ==0.) break;
+      hit =  (AliTPChit*) TPC->NextHit();
+    }
+    if (hit) {
+      hit0 = new AliTPChit(*hit); //Make copy of hit
+      hit = hit0;
+    }
+    else continue;
+    AliTPChit *hit1=(AliTPChit*)TPC->NextHit();       
+    if (hit1==0) continue;
+    if (hit1->fQ != 0.) continue;
+    Int_t i=hit->Track();
+    TParticle *p = (TParticle*)gAlice->Particle(i);
+    
+    if (p->GetFirstMother()>=0) continue;  //secondary particle
+    if (good[i] < good_number) continue;
+    if (p->Pt() < 0.1) continue;
+    if (TMath::Abs(p->Pz()/p->Pt())>0.999) continue;
+    printf("checking paricle %d with %d hits \n",i,good[i]);
+    
+    goodtracks[nt].eta = p->Eta();
+    goodtracks[nt].label=i;
+    goodtracks[nt].code=p->GetPdgCode();
+    goodtracks[nt].px=hit->X(); goodtracks[nt].py=hit->Y(); goodtracks[nt].pz=hit->Z();
+    goodtracks[nt].pt = p->Pt();
+    goodtracks[nt].nhits = good[i];
+    nt++;
+    
+    if (hit0) delete hit0;
+  }
+  
+  ofstream out(fname);
+  if(out) 
+    {
+      for (Int_t ngd=0; ngd<nt; ngd++)            
+       out<<goodtracks[ngd].label<<' '<<goodtracks[ngd].code<<' '<<
+         goodtracks[ngd].px<<' '<<goodtracks[ngd].py<<' '<<goodtracks[ngd].pz<<' '<<
+         goodtracks[ngd].pt<<' '<<goodtracks[ngd].eta<<' '<<goodtracks[ngd].nhits<<endl; 
+      
+    } 
+  else
+    printf("cannot open file\n");
+  
+  out.close();
+  
+  delete [] good;
+
+  evfile->Close();
+  digfile->Close();
+
+  delete transform;
+
+}
diff --git a/HLT/hough/GetGoodParticles.h b/HLT/hough/GetGoodParticles.h
new file mode 100644 (file)
index 0000000..f9d7bb7
--- /dev/null
@@ -0,0 +1,18 @@
+#ifndef GET_GOOD_PARTICLES
+#define GET_GOOD_PARTICLES
+
+struct GoodTrack 
+{
+
+  Int_t label;
+  Double_t eta;
+  Int_t code;
+  Double_t px,py,pz;
+  Double_t pt;
+  Int_t nhits;
+};
+
+
+void GetGoodParticles(Int_t,Char_t *,Char_t *);
+
+#endif