+void AliEMCALJetFinder::FillFromParticles()
+{
+// 26-feb-2002 PAI - for checking all chain
+// Work on particles level; accept all particle (not neutrino )
+
+ Double_t PX=0, PY=0, PZ=0, E=0; // checking conservation law
+ fNChTpc = 0;
+
+ ResetMap();
+ if (!fLego) BookLego();
+ fLego->Reset();
+//
+// Access particles information
+ Int_t npart = (gAlice->GetHeader())->GetNprimary();
+ if (fDebug >= 2 || npart<=0) {
+ printf(" AliEMCALJetFinder::FillFromParticles : npart %i\n", npart);
+ if(npart<=0) return;
+ }
+ fNt = npart;
+ fNtS = 0;
+ RearrangeParticlesMemory(npart);
+
+// Go through the particles
+ Int_t mpart, child1, child2, geantPdg;
+ Float_t pT, phi, eta, e=0, px=0, py=0, pz=0;
+ TParticle *MPart=0;
+ for (Int_t part = 0; part < npart; part++) {
+
+ fTrackList[part] = 0;
+
+ MPart = gAlice->Particle(part);
+ mpart = MPart->GetPdgCode();
+ child1 = MPart->GetFirstDaughter();
+ child2 = MPart->GetLastDaughter();
+ pT = MPart->Pt();
+ phi = MPart->Phi();
+ eta = MPart->Eta();
+
+ px = MPart->Px();
+ py = MPart->Py();
+ pz = MPart->Pz();
+ e = MPart->Energy();
+
+// see pyedit in Pythia's text
+ geantPdg = mpart;
+// Int_t kc = pycomp_(&mpart);
+// TString name = GetPythiaParticleName(mpart);
+ // printf(" mpart %6.6i;kc %6.6i -> gid %3.3i",mpart,kc,geantPdg);
+ //printf(" (%s)\n", name.Data());
+ if (IsThisPartonsOrDiQuark(mpart)) continue;
+ printf("%5i: %5i(%2i) px %5.1f py %5.1f pz %6.1f e %6.1f childs %5i,%5i \n",
+ part, mpart, geantPdg, px, py, pz, e, child1, child2);
+
+// exclude partons (21 - gluon, 92 - string)
+
+
+// exclude neutrinous also ??
+ if (fDebug >= 11 && pT>0.01)
+ printf("\n part:%5d mpart %5d eta %9.2f phi %9.2f pT %9.2f ",
+ part, mpart, eta, phi, pT);
+
+ fPtT[part] = pT;
+ fEtaT[part] = eta;
+ fPhiT[part] = phi;
+ fPdgT[part] = mpart;
+
+// final state only
+ if (child1 >= 0 && child1 < npart) continue;
+
+ // printf("%4i -> %5i(%3i) px %6.1f py %6.1f pz %7.1f e %8.2f child1 %5i %s\n",
+ // part, mpart, geantPdg, px, py, pz, e, child1, name.Data());
+ PX += px;
+ PY += py;
+ PZ += pz;
+ E += e;
+
+
+ if (TMath::Abs(eta) <= 0.9) fNChTpc++;
+// acceptance cut
+ if (TMath::Abs(eta) > 0.7) continue;
+ if (phi*180./TMath::Pi() > 120.) continue;
+//
+ if(fK0N==0 ) { // exclude neutral hadrons
+ if (mpart == kNeutron || mpart == kNeutronBar || mpart == kK0Long) continue;
+ }
+ fTrackList[part] = 1;
+ fLego->Fill(eta, phi, pT);
+
+ } // primary loop
+ printf("\n PX %8.2f PY %8.2f PZ %8.2f E %8.2f \n",
+ PX, PY, PZ, E);
+ DumpLego();
+ if(fhChPartMultInTpc) fhChPartMultInTpc->Fill(fNChTpc);
+}
+
+void AliEMCALJetFinder::FillFromPartons()
+{
+// function under construction - 13-feb-2002 PAI
+
+ if (fDebug >= 2)
+ printf("\n AliEMCALJetFinder::FillFromPartons()\n");
+ //
+
+ ResetMap();
+ if (!fLego) BookLego();
+ fLego->Reset();
+//
+// Access particle information
+ Int_t npart = (gAlice->GetHeader())->GetNprimary();
+ if (fDebug >= 2 || npart<=0)
+ printf("\n AliEMCALJetFinder::FillFromPartons : npart %i\n", npart);
+ fNt = 0; // for FindTracksInJetCone
+ fNtS = 0;
+
+// Go through the partons
+ Int_t statusCode=0;
+ for (Int_t part = 8; part < npart; part++) {
+ TParticle *MPart = gAlice->Particle(part);
+ Int_t mpart = MPart->GetPdgCode();
+ // Int_t child1 = MPart->GetFirstDaughter();
+ Float_t pT = MPart->Pt();
+ // Float_t p = MPart->P();
+ Float_t phi = MPart->Phi();
+ Float_t eta = MPart->Eta();
+ // Float_t theta = MPart->Theta();
+ statusCode = MPart->GetStatusCode();
+
+// accept partons (21 - gluon, 92 - string)
+ if (!(TMath::Abs(mpart) <= 6 || mpart == 21 ||mpart == 92)) continue;
+ if (fDebug > 1 && pT>0.01)
+ printf("\n part:%5d mpart %5d status %5d eta %8.2f phi %8.2f pT %8.2f ",
+ part, mpart, statusCode, eta, phi, pT);
+ // if (fDebug >= 3) MPart->Print();
+// accept partons before fragmentation - p.57 in Pythia manual
+// if(statusCode != 1) continue;
+// acceptance cut
+ if (TMath::Abs(eta) > 0.7) continue;
+ if (phi*180./TMath::Pi() > 120.) continue;
+// final state only
+// if (child1 >= 0 && child1 < npart) continue;
+//
+//
+ fLego->Fill(eta, phi, pT);
+
+ } // primary loop
+ DumpLego();
+}
+
+void AliEMCALJetFinder::BuildTrackFlagTable() {
+
+// Method to generate a lookup table for TreeK particles
+// which are linked to hits in the EMCAL
+//
+// --Author: J.L. Klay
+//
+// Access hit information
+ AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
+
+ TTree *TK = gAlice->TreeK(); // Get the number of entries in the kine tree
+ Int_t nKTrks = (Int_t) TK->GetEntries(); // (Number of particles created somewhere)
+
+ if(fTrackList) delete[] fTrackList; //Make sure we get rid of the old one
+ fTrackList = new Int_t[nKTrks]; //before generating a new one
+
+ for (Int_t i = 0; i < nKTrks; i++) { //Initialize members to 0
+ fTrackList[i] = 0;
+ }
+
+ TTree *treeH = gAlice->TreeH();
+ Int_t ntracks = (Int_t) treeH->GetEntries();
+//
+// Loop over tracks
+//
+ Int_t nbytes = 0;
+
+ for (Int_t track=0; track<ntracks;track++) {
+ gAlice->ResetHits();
+ nbytes += treeH->GetEvent(track);
+ if (pEMCAL) {
+
+//
+// Loop over hits
+//
+ for(AliEMCALHit* mHit=(AliEMCALHit*) pEMCAL->FirstHit(-1);
+ mHit;
+ mHit=(AliEMCALHit*) pEMCAL->NextHit())
+ {
+ Int_t iTrk = mHit->Track(); // track number
+ Int_t idprim = mHit->GetPrimary(); // primary particle
+
+ //Determine the origin point of this particle - it made a hit in the EMCAL
+ TParticle *trkPart = gAlice->Particle(iTrk);
+ TParticlePDG *trkPdg = trkPart->GetPDG();
+ Int_t trkCode = trkPart->GetPdgCode();
+ Double_t trkChg;
+ if (trkCode < 10000) { //Big Ions cause problems for
+ trkChg = trkPdg->Charge(); //this function. Since they aren't
+ } else { //likely to make it very far, set
+ trkChg = 0.0; //their charge to 0 for the Flag test
+ }
+ Float_t vxTrk = trkPart->Vx();
+ Float_t vyTrk = trkPart->Vy();
+ Float_t vrTrk = TMath::Sqrt(vxTrk*vxTrk+vyTrk*vyTrk);
+ fTrackList[iTrk] = SetTrackFlag(vrTrk,trkCode,trkChg);
+
+ //Loop through the ancestry of the EMCAL entrance particles
+ Int_t ancestor = trkPart->GetFirstMother(); //Get track's Mother
+ while (ancestor != -1) {
+ TParticle *ancPart = gAlice->Particle(ancestor); //get particle info on ancestor
+ TParticlePDG *ancPdg = ancPart->GetPDG();
+ Int_t ancCode = ancPart->GetPdgCode();
+ Double_t ancChg;
+ if (ancCode < 10000) {
+ ancChg = ancPdg->Charge();
+ } else {
+ ancChg = 0.0;
+ }
+ Float_t vxAnc = ancPart->Vx();
+ Float_t vyAnc = ancPart->Vy();
+ Float_t vrAnc = TMath::Sqrt(vxAnc*vxAnc+vyAnc*vyAnc);
+ fTrackList[ancestor] = SetTrackFlag(vrAnc,ancCode,ancChg);
+ ancestor = ancPart->GetFirstMother(); //Get the next ancestor
+ }
+
+ //Determine the origin point of the primary particle associated with the hit
+ TParticle *primPart = gAlice->Particle(idprim);
+ TParticlePDG *primPdg = primPart->GetPDG();
+ Int_t primCode = primPart->GetPdgCode();
+ Double_t primChg;
+ if (primCode < 10000) {
+ primChg = primPdg->Charge();
+ } else {
+ primChg = 0.0;
+ }
+ Float_t vxPrim = primPart->Vx();
+ Float_t vyPrim = primPart->Vy();
+ Float_t vrPrim = TMath::Sqrt(vxPrim*vxPrim+vyPrim*vyPrim);
+ fTrackList[idprim] = SetTrackFlag(vrPrim,primCode,primChg);
+ } // Hit Loop
+ } //if (pEMCAL)
+ } // Track Loop
+}
+
+Int_t AliEMCALJetFinder
+::SetTrackFlag(Float_t radius, Int_t code, Double_t charge) {
+
+ Int_t flag = 0;
+ Int_t parton = 0;
+ Int_t neutral = 0;
+
+ if (charge == 0) neutral = 1;
+
+ if (TMath::Abs(code) <= 6 ||
+ code == 21 ||
+ code == 92) parton = 1;
+
+ //It's not a parton, it's charged and it went through the TPC
+ if (!parton && !neutral && radius <= 84.0) flag = 1;
+
+ return flag;
+}
+
+
+
+void AliEMCALJetFinder::SaveBackgroundEvent()
+{
+// Saves the eta-phi lego and the tracklist
+//
+ if (fLegoB) {
+ fLegoB->Reset();
+ (*fLegoB) = (*fLegoB) + (*fLego);
+ if(fDebug)
+ printf("\n AliEMCALJetFinder::SaveBackgroundEvent() (fLegoB) %f = %f(fLego) \n",
+ fLegoB->Integral(), fLego->Integral());
+ }
+
+ if (fPtB) delete[] fPtB;
+ if (fEtaB) delete[] fEtaB;
+ if (fPhiB) delete[] fPhiB;
+ if (fPdgB) delete[] fPdgB;
+ if (fTrackListB) delete[] fTrackListB;
+
+ fPtB = new Float_t[fNtS];
+ fEtaB = new Float_t[fNtS];
+ fPhiB = new Float_t[fNtS];
+ fPdgB = new Int_t [fNtS];
+ fTrackListB = new Int_t [fNtS];
+
+ fNtB = 0;
+
+ for (Int_t i = 0; i < fNt; i++) {
+ if (!fTrackList[i]) continue;
+ fPtB [fNtB] = fPtT [i];
+ fEtaB[fNtB] = fEtaT[i];
+ fPhiB[fNtB] = fPhiT[i];
+ fPdgB[fNtB] = fPdgT[i];
+
+ fTrackListB[fNtB] = 1;
+ fNtB++;
+ }
+ fBackground = 1;
+ printf(" fNtB %i => fNtS %i #particles %i \n", fNtB, fNtS, fNt);
+}
+
+void AliEMCALJetFinder::InitFromBackground()
+{
+//
+//
+ if (fDebug) printf("\n AliEMCALJetFinder::InitFromBackground() ");
+
+ if (fLego) {
+ fLego->Reset();
+ (*fLego) = (*fLego) + (*fLegoB);
+ if(fDebug)
+ printf("\n AliEMCALJetFinder::SaveBackgroundEvent() (fLego) %f = %f(fLegoB) \n",
+ fLego->Integral(), fLegoB->Integral());
+ } else {
+ printf(" => fLego undefined \n");
+ }
+}
+
+