ROOT logo
/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to usec, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

//---------------------------------------------------------------------
// Jet Finder based on CDF algortihm
// Charged jet evolution and the underlying event in proton-antiproton collisions at 1.8 TeV
// Physical Review D, vol. 65, Issue 9, id. 092002
// http://www.phys.ufl.edu/~rfield/cdf/chgjet/chgjet_intro.html
// Authors : Adrian.Sevcenco@cern.ch (adriansev@spacescience.ro )
//           Daniel.Felea@cern.ch    (dfelea@spacescience.ro)
//           Ciprian.Mihai.Mitu@cern.ch  (mcm@spacescience.ro)
//           magali.estienne@subatech.in2p3.fr &
//           alexandre.shabetai@cern.ch (Modification of the input object (reader/finder splitting))
// ** 2011
// Modified accordingly to reader/finder splitting and new handling of neutral information
//---------------------------------------------------------------------

#include <Riostream.h>
#include <TMath.h>
#include <TBits.h>
#include <TFile.h>
#include <TH1F.h>
#include <TProfile.h>
#include <TVector2.h>

#include "AliAODJet.h"
#include "AliJetFinder.h"
#include "AliJetCalTrk.h"
#include "AliCdfJetFinder.h"
#include "AliCdfJetHeader.h"

ClassImp(AliCdfJetFinder)

//______________________________________________________________________________
AliCdfJetFinder::AliCdfJetFinder():
  AliJetFinder(),
  fHistos(0),
  fAODwrite(0),
  fAODtracksWrite(0),
  fAnalyseJets(0),
  fNJets(0),
  fNPart(0),
  fNInC(0),
  fNInN(0),
  fRadius(0.7),
  fMinJetParticles(1),
  fJetPtCut(0.),
  fVectParticle(NULL),
  fVectJet(NULL),
  fPtArray(NULL),
  fIdxArray(NULL)
{  
  // Default constructor
}

//______________________________________________________________________________
AliCdfJetFinder::~AliCdfJetFinder()
{
  // destructor
  Clean();

}

//______________________________________________________________________________
void AliCdfJetFinder::CreateOutputObjects(TList * const histos)
{
  // Create the list of histograms. Only the list is owned.
  fHistos = histos;

  //  gStyle->SetOptStat(11111111);

  TH1F *h1 = new TH1F ("histo1", "Pt distribution of jets", 200, 0,200);  // 1GeV/bin
  h1->SetStats(kTRUE);
  h1->GetXaxis()->SetTitle("P_{T} of jets");
  h1->GetYaxis()->SetTitle("Number of jets");
  h1->GetXaxis()->SetTitleColor(1);
  h1->SetMarkerStyle(kFullCircle);
  fHistos->Add(h1);

  TH1F *h2 = new TH1F ("histo2", "Eta distribution of jets", 240, -1.2,1.2); // 1 unit of rapidity / 100 bin
  h2->SetStats(kTRUE);
  h2->GetXaxis()->SetTitle("Eta of jets");
  h2->GetYaxis()->SetTitle("Number of jets");
  h2->GetXaxis()->SetTitleColor(1);
  h2->SetMarkerStyle(kFullCircle);
  fHistos->Add(h2);

  TH1F *h3 = new TH1F ("histo3", "Phi distribution of jets", 400, -4,4);
  h3->SetStats(kTRUE);
  h3->GetXaxis()->SetTitle("Phi of jets");
  h3->GetYaxis()->SetTitle("Number of jets");
  h3->GetXaxis()->SetTitleColor(1);
  h3->SetMarkerStyle(kFullCircle);
  fHistos->Add(h3);

  TH1F *h4 = new TH1F ("histo4", "Multiplicity of jets", 40, 0,40);  // 1 unit of multiplicity /bin
  h4->SetStats(kTRUE);
  h4->GetXaxis()->SetTitle("Particles in jets");
  h4->GetYaxis()->SetTitle("Number of jets");
  h4->GetXaxis()->SetTitleColor(1);
  h4->SetMarkerStyle(kFullCircle);
  fHistos->Add(h4);

  TH1F *h5 = new TH1F ("histo5", "Distribution of jets in events", 100, 0,100);
  h5->SetStats(kTRUE);
  h5->GetXaxis()->SetTitle("Number of jets");
  h5->GetYaxis()->SetTitle("Number of events");
  h5->GetXaxis()->SetTitleColor(1);
  h5->SetMarkerStyle(kFullCircle);
  fHistos->Add(h5);

  TH1F *h6 = new TH1F ("histo6", "Jet1 Charged Multiplicity Distribution", 30, 0,30);
  h6->SetStats(kTRUE);
  h6->GetXaxis()->SetTitle("N_{chg}");
  h6->GetYaxis()->SetTitle("Number of jets");
  h6->GetXaxis()->SetTitleColor(1);
  h6->SetMarkerStyle(kFullCircle);
  fHistos->Add(h6);

  TProfile * h7 = new TProfile ("histo7","N_{chg}(jet1) vs P_{T}(charged jet1)", 200, 0. ,200. , 0.,200. ) ;
  h7->SetStats(kTRUE);
  h7->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h7->GetYaxis()->SetTitle("<N_{chg}(jet1)> in 1 GeV/c bin");
  h7->GetXaxis()->SetTitleColor(1);
  h7->SetMarkerStyle(kFullCircle);
  fHistos->Add(h7);

  TH1F *h8 = new TH1F ("histo8", "Charge momentum distribution for leading jet", 120, 0 , 1.2);
  h8->SetStats(kTRUE);
  h8->GetXaxis()->SetTitle("Jets");
  h8->GetYaxis()->SetTitle("Particle distribution");
  h8->GetXaxis()->SetTitleColor(1);
  h8->SetMarkerStyle(kFullCircle);
  fHistos->Add(h8);

  TProfile *h9 = new TProfile ("histo9", "N_{chg} vs the Azimuthal Angle from Charged Jet1", 50 , 0. , 180. , 0 , 20 );
  h9->SetStats(kTRUE);
  h9->GetXaxis()->SetTitle("#Delta#phi (degrees)");
  h9->GetYaxis()->SetTitle("<N_{chg}> in 3.6 degree bin");
  h9->GetXaxis()->SetTitleColor(1);
  h9->SetMarkerStyle(kFullCircle);
  fHistos->Add(h9);

  TProfile *h10 = new TProfile ("histo10", "P_{T} sum vs the Azimuthal Angle from Charged Jet1", 50 , 0. , 180. , 0 , 100 );
  h10->SetStats(kTRUE);
  h10->GetXaxis()->SetTitle("#Delta#phi (degrees)");
  h10->GetYaxis()->SetTitle("<P_{T} sum> in 3.6 degree bin");
  h10->GetXaxis()->SetTitleColor(1);
  h10->SetMarkerStyle(kFullCircle);
  fHistos->Add(h10);

  TH1F *h11 = new TH1F ("histo11", " \"Transverse\" Pt Distribution ", 70, 0 , 14);
  h11->SetStats(kTRUE);
  h11->GetXaxis()->SetTitle("P_{T} (GeV/c)");
  h11->GetYaxis()->SetTitle("dN_{chg}/dP_{T} (1/GeV/c)");
  h11->GetXaxis()->SetTitleColor(1);
  h11->SetMarkerStyle(kFullCircle);
  fHistos->Add(h11);

  TH1F *h20 = new TH1F ("histo20", "Distribution of R in leading jet", 400, 0.,4.);
  h20->SetStats(kTRUE);
  h20->GetXaxis()->SetTitle("R [formula]");
  h20->GetYaxis()->SetTitle("dN/dR");
  h20->GetXaxis()->SetTitleColor(1);
  h20->SetMarkerStyle(kFullCircle);
  fHistos->Add(h20);

  TProfile * h21 = new TProfile ("histo21","N_{chg}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0., 50. , 0., 30. ) ;
  h21->SetStats(kTRUE);
  h21->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h21->GetYaxis()->SetTitle("<N_{chg}(in the event - including jet1)> in 1 GeV/c bin");
  h21->GetXaxis()->SetTitleColor(1);
  h21->SetMarkerStyle(kFullCircle);
  fHistos->Add(h21);

  TProfile * h22 = new TProfile ("histo22","PT_{sum}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0. , 50. , 0., 50. ) ;
  h22->SetStats(kTRUE);
  h22->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h22->GetYaxis()->SetTitle("<PT_{sum}(in the event - including jet1)> in 1 GeV/c bin");
  h22->GetXaxis()->SetTitleColor(1);
  h22->SetMarkerStyle(kFullCircle);
  fHistos->Add(h22);

  TProfile * h21Toward = new TProfile ("histo21_toward","N_{chg}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0., 50. , 0., 12. ) ;
  h21Toward->SetStats(kTRUE);
  h21Toward->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h21Toward->GetYaxis()->SetTitle("<N_{chg}(in the event - including jet1)> in 1 GeV/c bin");
  h21Toward->GetXaxis()->SetTitleColor(1);
  h21Toward->SetMarkerStyle(kFullCircle);
  fHistos->Add(h21Toward);

  TProfile * h21Transverse = new TProfile ("histo21_transverse","N_{chg}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0., 50. , 0., 12. ) ;
  h21Transverse->SetStats(kTRUE);
  h21Transverse->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h21Transverse->GetYaxis()->SetTitle("<N_{chg}(in the event - including jet1)> in 1 GeV/c bin");
  h21Transverse->GetXaxis()->SetTitleColor(1);
  h21Transverse->SetMarkerStyle(kFullCircle);
  fHistos->Add(h21Transverse);

  TProfile * h21Away = new TProfile ("histo21_away","N_{chg}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0., 50. , 0., 12. ) ;
  h21Away->SetStats(kTRUE);
  h21Away->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h21Away->GetYaxis()->SetTitle("<N_{chg}(in the event - including jet1)> in 1 GeV/c bin");
  h21Away->GetXaxis()->SetTitleColor(1);
  h21Away->SetMarkerStyle(kFullCircle);
  fHistos->Add(h21Away);

  TProfile * h22Toward = new TProfile ("histo22_toward","PT_{sum}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0. , 50. , 0., 50. ) ;
  h22Toward->SetStats(kTRUE);
  h22Toward->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h22Toward->GetYaxis()->SetTitle("<PT_{sum}(in the event - including jet1)> in 1 GeV/c bin");
  h22Toward->GetXaxis()->SetTitleColor(1);
  h22Toward->SetMarkerStyle(kFullCircle);
  fHistos->Add(h22Toward);

  TProfile * h22Transverse = new TProfile ("histo22_transverse","PT_{sum}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0. , 50. , 0., 50. ) ;
  h22Transverse->SetStats(kTRUE);
  h22Transverse->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h22Transverse->GetYaxis()->SetTitle("<PT_{sum}(in the event - including jet1)> in 1 GeV/c bin");
  h22Transverse->GetXaxis()->SetTitleColor(1);
  h22Transverse->SetMarkerStyle(kFullCircle);
  fHistos->Add(h22Transverse);

  TProfile * h22Away = new TProfile ("histo22_away","PT_{sum}(in the event - including jet1) vs P_{T}(charged jet1)", 200, 0. , 50. , 0., 50. ) ;
  h22Away->SetStats(kTRUE);
  h22Away->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h22Away->GetYaxis()->SetTitle("<PT_{sum}(in the event - including jet1)> in 1 GeV/c bin");
  h22Away->GetXaxis()->SetTitleColor(1);
  h22Away->SetMarkerStyle(kFullCircle);
  fHistos->Add(h22Away);

  TH1F *h23Toward = new TH1F ("histo23_toward","'Toward' Pt Distribution of charged particles", 200, 0., 14.);
  h23Toward->SetStats(kTRUE);
  h23Toward->GetXaxis()->SetTitle("P_{T} (charged) (GeV/c)");
  h23Toward->GetYaxis()->SetTitle("dN_{chg}/dP_{T} (1/GeV/c)");
  h23Toward->GetXaxis()->SetTitleColor(1);
  h23Toward->SetMarkerStyle(kFullCircle);
  fHistos->Add(h23Toward);

  TH1F *h23Transverse = new TH1F ("histo23_transverse","'Transverse' Pt Distribution of charged particles", 200, 0., 14.);
  h23Transverse->SetStats(kTRUE);
  h23Transverse->GetXaxis()->SetTitle("P_{T} (charged) (GeV/c)");
  h23Transverse->GetYaxis()->SetTitle("dN_{chg}/dP_{T} (1/GeV/c)");
  h23Transverse->GetXaxis()->SetTitleColor(1);
  h23Transverse->SetMarkerStyle(kFullCircle);
  fHistos->Add(h23Transverse);

  TH1F *h23Away = new TH1F ("histo23_away","'Away' Pt Distribution of charged particles", 200, 0., 14.);
  h23Away->SetStats(kTRUE);
  h23Away->GetXaxis()->SetTitle("P_{T} (charged) (GeV/c)");
  h23Away->GetYaxis()->SetTitle("dN_{chg}/dP_{T} (1/GeV/c)");
  h23Away->GetXaxis()->SetTitleColor(1);
  h23Away->SetMarkerStyle(kFullCircle);
  fHistos->Add(h23Away);

  TProfile * h24 = new TProfile ("histo24","Jet1 Size vs P_{T}(charged jet1)", 200, 0., 50. , 0., 0.5) ;
  h24->SetStats(kTRUE);
  h24->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h24->GetYaxis()->SetTitle("<R(chgjet1)> in 1 GeV/c bin");
  h24->GetXaxis()->SetTitleColor(1);
  h24->SetMarkerStyle(kFullCircle);
  fHistos->Add(h24);

  TProfile * h25 = new TProfile ("histo25","Jet1 Size vs P_{T}(charged jet1)", 200, 0., 50. , 0., 0.5) ;
  h25->SetStats(kTRUE);
  h25->GetXaxis()->SetTitle("P_{T} (charged jet1) (GeV/c)");
  h25->GetYaxis()->SetTitle("<R(chgjet1)> in 1 GeV/c bin");
  h25->GetXaxis()->SetTitleColor(1);
  h25->SetMarkerStyle(kFullCircle);
  fHistos->Add(h25);

  TProfile *h26 = new TProfile ("histo26", "N_{chg} vs the Distance R from Charged Jet1", 30, 0., 0.6, 0., 0.8);
  h26->SetStats(kTRUE);
  h26->GetXaxis()->SetTitle("Distance R");
  h26->GetYaxis()->SetTitle("<N_{chg}> in 0.02 bin");
  h26->GetXaxis()->SetTitleColor(1);
  h26->SetMarkerStyle(kFullCircle);
  fHistos->Add(h26);

  TProfile *h27 = new TProfile ("histo27", "N_{chg} vs the Distance R from Charged Jet1", 30, 0., 0.6, 0., 0.8);
  h27->SetStats(kTRUE);
  h27->GetXaxis()->SetTitle("Distance R");
  h27->GetYaxis()->SetTitle("<N_{chg}> in 0.02 bin");
  h27->GetXaxis()->SetTitleColor(1);
  h27->SetMarkerStyle(kFullCircle);
  fHistos->Add(h27);

  TProfile *h28 = new TProfile ("histo28", "PT_{sum} vs the Distance R from Charged Jet1", 30, 0., 0.6, 0.01, 10.);
  h28->SetStats(kTRUE);
  h28->GetXaxis()->SetTitle("Distance R");
  h28->GetYaxis()->SetTitle("<PT_{sum} (GeV/c)> in 0.02 bin");
  h28->GetXaxis()->SetTitleColor(1);
  h28->SetMarkerStyle(kFullCircle);
  fHistos->Add(h28);

  TProfile *h29 = new TProfile ("histo29", "PT_{sum} vs the Distance R from Charged Jet1", 30, 0., 0.6, 0.01, 10.);
  h29->SetStats(kTRUE);
  h29->GetXaxis()->SetTitle("Distance R");
  h29->GetYaxis()->SetTitle("<PT_{sum} (GeV/c)> in 0.02 bin");
  h29->GetXaxis()->SetTitleColor(1);
  h29->SetMarkerStyle(kFullCircle);
  fHistos->Add(h29);

}

//______________________________________________________________________________
void AliCdfJetFinder::FindJets()
{
  // Jet Algorithm:
  //  * Order all charged particles according to their PT.
  //  * Start with the highest PT particle and include in the "jet" all particles within the "radius" R = 0.7
  //     (considering each particle in the order of decreasing PT and recalculating the centroid of the jet after
  //     each new particle is added to the jet).
  //  * Go to the next highest PT particle (not already included in a jet) and include in the "jet" all particles
  //     (not already included in a jet) within the radius R =0.7.
  //  * Continue until all particles are in a "jet".
  if (fDebug) { printf("AliCDJetfinder::FindJets() %d \n", __LINE__ ); }
  AliCdfJetHeader *header = (AliCdfJetHeader*)fHeader;

  if (header)
    {
      fDebug            = header->GetDebug();
      fAODwrite         = header->IsAODwrite() ;       // write jets to AOD
      fAODtracksWrite   = header->IsAODtracksWrite() ; // write jet tracks to AOD
      fRadius           = header->GetRadius();      // get Radius from jet finder header
      fMinJetParticles  = header->GetMinPartJet (); // get minimum multiplicity of an jet
      fJetPtCut         = header->GetJetPtCut ();   // get minimum of jet pt
      fAnalyseJets      = header->GetAnalyseJets(); // get analyse jet
    }
  else
    { cout << "Header not found" << endl; return; }

  InitData();

  if (!fNPart) { 
    if (fDebug) {
      cout << "entries = 0 ; Event empty !!!" << endl ;
    }
    // no need to call clean, InitData does not 
    // create pointers if npart == 0
    return; 
  } // if event empty then exit

  FindCones();
 
  ComputeConesWeight();
 
  if (fAODwrite) { 
    if(fDebug)cout << "Writing AOD" << endl ; 
    WriteJets();
  }
 
  if (fAnalyseJets) AnalizeJets();
 
  Clean();

}

//______________________________________________________________________________
void AliCdfJetFinder::InitData()
{
  // initialisation of variables and data members

  if( fHeader->GetDebug() && fCalTrkEvent->GetNCalTrkTracks()  == 0) { cout << "No charged tracks found" << endl; }

  fNPart = fCalTrkEvent->GetNCalTrkTracks() ;

  if ( fCalTrkEvent->GetNCalTrkTracks() ) { return; } // if event empty then exit

  fVectParticle = new varContainer* [fNPart]; // container for Particles

  fPtArray  = new Double_t [fCalTrkEvent->GetNCalTrkTracks()] ; 
  fIdxArray = new Int_t    [fCalTrkEvent->GetNCalTrkTracks()] ; // index array of sorted pts

  // initialisation of momentum and index arrays
  for (  Int_t i = 0 ; i <fCalTrkEvent->GetNCalTrkTracks() ; i++ )
    {// SORTING STEP :: fPtArray with data from CalTrkTracks

      // INITIALISATION of local arrays for temporary storage
      varContainer *aParticle = new varContainer;
      aParticle->pt   = fCalTrkEvent->GetCalTrkTrack(i)->GetPt();
      aParticle->eta  = fCalTrkEvent->GetCalTrkTrack(i)->GetEta();
      aParticle->phi  = TVector2::Phi_mpi_pi ( fCalTrkEvent->GetCalTrkTrack(i)->GetPhi() ); // normalize to -pi,pi
      aParticle->njet = -999;

      fVectParticle[i] = aParticle;  // vector of Particles

      // initializing arrays
      fIdxArray [i] = -999 ;
      fPtArray [i] = aParticle->pt ;
    }

  TMath::Sort ( fNPart, fPtArray, fIdxArray ) ; // get a sorted array of indexes

}

//______________________________________________________________________________
void AliCdfJetFinder::FindCones()
{
  // parsing of particles in event and estlabish jets (label them with jet index)

  Double_t  ptSeed = 0. , etaSeed = 0. , phiSeed = 0. ; // leading particle params
  Double_t pttmp = 0. , etatmp = 0. , phitmp = 0. ; // temporary variables to be used in various calculations
  Double_t deta = 0. , dphi = 0. , dcomputed = 0. ;
  Bool_t injet = 0 ;

  fNJets = -1 ; // n jets in this event
  Int_t idxPtSort = -1 ;  // index of array of sorted pt indexes

  if (fDebug) { cout << "\n\n\n\n\n\n------------------\nBegin Event Analysis\n------------------\n\n" << endl ;}

  if(fDebug)cout << "fNPart = " << fNPart << endl;

  TBits lkupTable ( fNPart ) ;  // bit container ; 1-to-1 corespondence with fIdxArray

  while ( lkupTable.CountBits() != (UInt_t)fNPart )
    { // loop over particles in event until all flags are set
      UInt_t firstnonflagged = lkupTable.FirstNullBit() ; // set the index to the first NON flagged bit ; less conditions

      if(fDebug)cout << "\n\nfirst_non_flagged : " << firstnonflagged << endl;

      ++fNJets; // incrementing the jet counter
      if (fDebug) { printf("JET %d \n", fNJets); }

      ptSeed = 0. ; etaSeed = 0. ; phiSeed = 0. ;  // reseting leading particle params

      for (  UInt_t ipart = firstnonflagged ; ipart < (UInt_t)fNPart ; ipart++ )
	{// iteration over particles in event
	  // the loop is done over sorted array of pt
	  idxPtSort = fIdxArray[ipart] ;  // index of particle ! fIdxArray is an index list pt sorted

	  if ( lkupTable.TestBitNumber(ipart) ) { continue; } // if 4vector is already flagged skip it

	  //init computed and used vars
	  pttmp = 0. ; etatmp = 0. ; phitmp = 0. ;
	  deta = 0. ; dphi = 0. ; dcomputed = 0. ; injet = 0 ;

	  //taking info from fVectParticle ;
	  pttmp = fVectParticle[idxPtSort]->pt ;
	  etatmp = fVectParticle[idxPtSort]->eta ;
	  phitmp = fVectParticle[idxPtSort]->phi ;

	  if ( ipart == firstnonflagged )
	    {// this is first particle in event; leading particle
	      // begin the search around this particle in a fRadius

	      // CENTRE OF THE JET
	      ptSeed = pttmp ; etaSeed = etatmp ; phiSeed = phitmp ; // seeding the jet with first particle idxPtSort

	      lkupTable.SetBitNumber ( ipart ) ; // flag the index of particle in lkup_table
	      fVectParticle[idxPtSort]->njet = fNJets ; // associate particle with current jet number

	      if (fDebug) { printf("\nLeading particle :: particle index = %d ;  at sorted index %d ; in jet %d \n", idxPtSort, ipart, fNJets); }
	      if (fDebug) { printf("pt= %g ; eta= %g ; phi = %g \n", pttmp, etatmp, phitmp) ; }
	      if (fDebug) { lkupTable.Print() ;}

	      continue ; // skip to next particle
	    }

	  // condition to be in jet
	  deta = etatmp - etaSeed ;
	  dphi = TVector2::Phi_mpi_pi (phitmp - phiSeed) ; // computing dphi and normalizing to (0,2pi) interval in one step

	  dcomputed = TMath::Hypot(deta, dphi) ; // Distance(fRadius) to (eta,phi) seed

	  injet = ( ( fRadius - dcomputed ) >= 0.000000001 ) ? 1 : 0 ; // if r_computed is within jet_r in_jet == 1 else 0

	  if ( injet )
	    { // calculus of jet variables
	      lkupTable.SetBitNumber ( ipart ) ;  // flag the index of particle in lkup_table
	      fVectParticle[idxPtSort]->njet = fNJets ; // setting in particle list the associated jet

	      if (fDebug) { printf("\njet particle :: particle index = %d ; at sorted index %d ; in jet %d ; found at radius %g ;  \n", idxPtSort, ipart, fNJets, dcomputed); }
	      if (fDebug) { printf("pt= %g ; eta= %g ; phi = %g \n", pttmp, etatmp, phitmp) ; }
	      if (fDebug) { lkupTable.Print() ;}

	      continue ; // skip to next particle
	    }

	}
      // end of iteration over event; one jet definition of content ; jet parameters to be computed later
    }

}

//______________________________________________________________________________
void AliCdfJetFinder::ComputeConesWeight()
{
  // computing of jets Pt, Eta and Phi (centre of weight in (eta,phi) plane)
  // rescan the vector of particles by identify them by asociate jet number for computing of weight centre

  // JET CONTAINER
  fVectJet      = new varContainer* [fNJets]; // container for Jets

  Double_t ptJet, ptJet2 , etaJet , phiJet ; Int_t npartJet ;
  Double_t pttmp = 0. , etatmp = 0. , phitmp = 0. ; // temporary variables to be used in various calculations
  Int_t idxPtSort = -999 ;  // index of array of sorted pt indexes

  for(  Int_t jet = 0 ; jet < fNJets ; jet++ )
    {
      if (fDebug) { printf("\n\n--- Computing weight of Jet %d \n", jet ); }
      npartJet = 0 ; ptJet = 0. ; etaJet = 0. ; phiJet = 0. ; // reset variables for a new computation

      for (  Int_t ipart = 0 ; ipart < fNPart ; ipart++ )
	{// iteration over particles in event
	  // the loop is done over sorted array of pt
	  idxPtSort = fIdxArray[ipart] ;  // index of particle ! fIdxArray is an index list pt sorted

	  if ( fVectParticle[idxPtSort]->njet == jet )
	    {
	      ++npartJet; // incrementing the counter of jet particles

	      //taking info from fVectParticle ;
	      pttmp = fVectParticle[idxPtSort]->pt ;
	      etatmp = fVectParticle[idxPtSort]->eta ;
	      phitmp = TVector2::Phi_mpi_pi (fVectParticle[idxPtSort]->phi) ;

	      //      jet_new_angular_coordinate = jet_old_angular_coordinate * jet_old_pt / jet_new_pt +
	      //                                    part[i]_angular_coordinate * part[i]_pt/jet_new_pt

	      ptJet2 = ptJet + pttmp ;

	      etaJet = etaJet * ptJet / ptJet2 +  etatmp * pttmp / ptJet2 ;
	      phiJet = phiJet * ptJet / ptJet2 +  phitmp * pttmp / ptJet2 ;

	      ptJet = ptJet2 ;

	    }
	  // add a particle and recalculation of centroid
	}
      // end of 1 jet computation

      varContainer *aJet = new varContainer;  // Jet container
      aJet->pt = ptJet; aJet->eta = etaJet; aJet->phi = phiJet; aJet->njet = npartJet; // setting jet vars in container
      fVectJet[jet] = aJet;   // store the number of the jet(fNJets) and increment afterwards

      if (fDebug) { printf ("=== current jet %d : npartjet= %d ; pt_jet= %g ; eta_jet = %g ; phi_jet = %g \n\n\n",
			    jet,     npartJet,      ptJet,      etaJet,       phiJet ) ; }

    }
  //end loop over jets

}

//______________________________________________________________________________
void AliCdfJetFinder::WriteJets()  
{ 
  // Writing AOD jets and AOD tracks

  for(  Int_t jetnr = 0 ; jetnr < fNJets ; jetnr++ )
    {
      Double_t pt = 0., eta = 0., phi = 0., // jet variables
	px = 0., py = 0., pz = 0., en = 0.; // convert to 4-vector
      pt  = fVectJet[ jetnr ]->pt   ; // pt  of jet
      eta = fVectJet[ jetnr ]->eta  ; // eta of jet
      phi = fVectJet[ jetnr ]->phi  ; // phi of jet

      px = pt * TMath::Cos ( phi ) ;
      py = pt * TMath::Sin ( phi ) ;
      pz = pt / TMath::Tan ( 2.0 * TMath::ATan ( TMath::Exp ( -eta ) ) ) ;
      en = TMath::Sqrt ( px * px + py * py + pz * pz );

      AliAODJet jet (px, py, pz, en);


      if (fDebug) jet.Print("");

      if (fAODtracksWrite)
	{
	  for (  Int_t jetTrack = 0; jetTrack < fCalTrkEvent->GetNCalTrkTracks(); jetTrack++ )
	    {
	      // The first if condition below has to be checked
	      if ( fVectParticle[jetTrack]->njet == jetnr ) { jet.AddTrack(fCalTrkEvent->GetCalTrkTrack(jetTrack)->GetTrackObject()) ; }
	    }
	}
      // tracks REFs written in AOD
      AddJet(jet);

    }
  //jets vector parsed and written to AOD
}

//______________________________________________________________________________
void AliCdfJetFinder::AnalizeJets()
{
  // analyzing of jets and filling of histograms

  const Double_t kPI = TMath::Pi();
    
  //persistent pointer to histo20
  TH1F *hR = (TH1F*)fHistos->FindObject("histo20");

  Int_t   *jetsptidx = 0;     // sorted array of jets pt
  Double_t    *jetspt = 0;     // array of jets pts
  Int_t leadingjetindex = -1 ;   // index of leading jet from fVectJet
  Int_t partleadjet = 0 ; // number of particles in leading jet
  Double_t   ptleadjet = 0. ; // pt  of leading jet
  Double_t  etaleadjet = 0. ; // eta of leading jet
  Double_t  phileadjet = 0. ; // phi of leading jet

  jetsptidx = new Int_t    [fNJets] ;
  jetspt    = new Double_t [fNJets] ;

  //________________________________________________________________________________________
  //  Jet sorting and finding the leading jet that coresponds to cuts in pt and multiplicity
  //________________________________________________________________________________________

  // filing the idx_ptjets array
  if (fDebug) printf("List of unsorted jets:\n");
  for(  Int_t i = 0 ; i < fNJets ; i++ )
    {
      jetsptidx [i] = 0 ;
      jetspt    [i] = fVectJet[i]->pt ;
      if (fDebug) { cout << "   jet found: " << i << " npartjet=" << fVectJet[i]->njet << " ; jets_pt = " << jetspt[i] << endl; }
    }

  TMath::Sort ( fNJets, jetspt , jetsptidx ) ; // sorting pt of jets

  // selection of leading jet
  // looping over jets searching for __first__ one that coresponds to cuts
  for(  Int_t i = 0 ; i < fNJets ; i++ )
    {
      if ( ( fVectJet[ jetsptidx[i] ]->njet >= fMinJetParticles ) && ( fVectJet[ jetsptidx[i] ]->pt >= fJetPtCut ) )
	{
	  leadingjetindex = jetsptidx[i] ;
	  partleadjet = fVectJet[ leadingjetindex ]->njet ; // number of particles in leading jet
	  ptleadjet = fVectJet[ leadingjetindex ]->pt   ; // pt  of leading jet
	  etaleadjet = fVectJet[ leadingjetindex ]->eta  ; // eta of leading jet
	  phileadjet = fVectJet[ leadingjetindex ]->phi  ; // phi of leading jet

	  if (fDebug)
	    { printf("Leading jet %d : npart= %d ; pt= %g ; eta = %g ; phi = %g \n", leadingjetindex, partleadjet, ptleadjet, etaleadjet, phileadjet ); }

	  break ;
	}
    }
  // end of selection of leading jet



  //////////////////////////////////////////////////
  ////  Computing of values used in histograms
  //////////////////////////////////////////////////

  //___________________________________________________________________________
  // pt_sum of all particles in event
  //___________________________________________________________________________
  if (fDebug) cout << "Computing sum of pt in event" << endl ;
  Double_t ptsumevent = 0.;
  for (  Int_t i = 0 ; i< fNPart ; i++ ) { ptsumevent += fVectParticle[i]->pt ; }
  if (fDebug) printf ("Sum of all Pt in event : pt_sum_event = %g", ptsumevent) ;

  //___________________________________________________________________________
  // Filling an array with indexes of leading jet particles
  //___________________________________________________________________________
  Int_t * idxpartLJ = new Int_t [partleadjet] ;
  Int_t counterpartleadjet = 0;

  if (fDebug) cout << "Filling an array with indexes of leading jet particles" << endl;

  for( Int_t i = 0 ; i < fNPart ; i++ )
    {
      if ( fVectParticle[i]->njet == leadingjetindex )
	{  idxpartLJ[counterpartleadjet++] = i ; }
    }

  if ( (counterpartleadjet-1) > partleadjet ) { cout << " Counter_part_lead_jet > part_leadjet !!!!" << endl;}


  //___________________________________________________________________________
  // Calculus of part distribution in leading jet
  //___________________________________________________________________________
  Double_t z = 0. ;
  Double_t *zpartljet = new Double_t [ partleadjet ] ; // array of z of particles in leading jet

  if (fDebug) cout << "Entering loop of calculus of part distribution in leading jet" << endl ;

  for( Int_t j = 0 ; j < partleadjet ; j++ )
    {
      Double_t zj = fVectParticle[idxpartLJ[j]]->pt ;
      z =  zj / ptleadjet ;
      zpartljet [j] = z ;
      if (fDebug) cout << "idx_leadjet_part[j] = " << idxpartLJ[j]
		       << " p of particle = " << zj
		       << " pt lead jet = " << ptleadjet
		       << " Z = " << z << endl;
    }


  //___________________________________________________________________________
  // array of delta phi's between phi of particles and leading jet phi
  //___________________________________________________________________________
  if (fDebug) cout << "array of delta phi's between phi of particles and leading jet phi" << endl;
  Double_t dphipartLJ = 0. ;
  Double_t *dphipartljet = new Double_t [fNPart];
  for(  Int_t part = 0 ; part < fNPart ; part++ )
    {
      dphipartLJ = fVectParticle[part]->phi - phileadjet ;
      dphipartLJ = TVector2::Phi_mpi_pi (dphipartLJ) ; // restrict the delta phi to (-pi,pi) interval
      dphipartljet [part] = dphipartLJ ;
      if (fDebug) printf("part= %d ; dphi_partLJ = %g  \n", part, dphipartLJ );
    }


  //______________________________________________________________________________
  //  Pt distribution for all particles
  //______________________________________________________________________________
  TH1F * hpt = (TH1F*)fHistos->FindObject("histo11");
  if ( hpt ) { for (  Int_t i = 0 ; i < fNPart ; i++ ) { hpt->Fill( fVectParticle[i]->pt ); } }

  //___________________________________________________________________________
  // Recomputing of radius of particles in leading jet
  //___________________________________________________________________________
  if (fDebug) { printf("   Searching particles with jet index %d\n", leadingjetindex); }

  Double_t ddeta = 0. , ddphi = 0. , rpart = 0. ;

  for( Int_t j = 0 ; j < partleadjet ; j++ )
    {
      ddeta = etaleadjet - fVectParticle[idxpartLJ[j]]->eta;

      Double_t phitmp = fVectParticle[idxpartLJ[j]]->phi ;

      ddphi = TVector2::Phi_mpi_pi ( phileadjet - phitmp ) ; // restrict the delta phi to (-pi,pi) interval

      rpart = TMath::Hypot (ddeta, ddphi) ;

      if (fDebug) printf ("Particle %d with Re-Computed radius = %f ", idxpartLJ[j], rpart) ;
      if ( (rpart - fRadius) >= 0.00000001 )
	{ if (fDebug) printf ("    bigger than selected radius of %f\n", fRadius ); }
      else
	{ if (fDebug) printf ("\n") ; }

      if (hR) hR->Fill(rpart);

    }



  //_______________________________________________________________________
  // Computing of radius that contain 80% of Leading Jet ( PT and multiplicity )
  //_______________________________________________________________________
  Double_t corepartleadjet = 0.8 * partleadjet ;
  Double_t coreptleadjet = 0.8 * ptleadjet ;
  Int_t countercorepart = 0 ;
  Double_t countercorept = 0. ;
  Int_t sortedindex = -1 ;

  TProfile * hprof24 = (TProfile*)fHistos->FindObject("histo24");
  TProfile * hprof25 = (TProfile*)fHistos->FindObject("histo25");
  TProfile * hprof26 = (TProfile*)fHistos->FindObject("histo26");
  TProfile * hprof27 = (TProfile*)fHistos->FindObject("histo27");
  TProfile * hprof28 = (TProfile*)fHistos->FindObject("histo28");
  TProfile * hprof29 = (TProfile*)fHistos->FindObject("histo29");


  if ((hprof24) && (hprof25) && (hprof26) && (hprof27) && (hprof28) && (hprof29) )
    {
      for(  Int_t part = 0 ; part < fNPart ; part++ )
	{
	  Double_t pttmp = 0. ; Double_t etatmp = 0. ; Double_t phitmp = 0. ; // temporary variables
	  Double_t dpart = 0. ;
	  sortedindex = fIdxArray[part] ;

	  if ( fVectParticle [ sortedindex ]->njet == leadingjetindex )
	    {
	      pttmp = fVectParticle[sortedindex]->pt ;
	      etatmp = fVectParticle[sortedindex]->eta ;
	      phitmp = fVectParticle[sortedindex]->phi ;

	      ++countercorepart ;
	      countercorept += pttmp ;

	      dpart = TMath::Hypot ( etaleadjet - etatmp, TVector2::Phi_mpi_pi (phileadjet - phitmp) ) ;

	      if ( countercorepart <=  corepartleadjet ) { hprof24->Fill(ptleadjet, dpart); }
	      if ( countercorept <= coreptleadjet ) { hprof25->Fill(ptleadjet, dpart); }

	      if (ptleadjet >  5.) { hprof26->Fill(dpart, countercorepart); hprof28->Fill(dpart, countercorept); }
	      if (ptleadjet > 30.) { hprof27->Fill(dpart, countercorepart); hprof29->Fill(dpart, countercorept); }

	    }
	}
    }

  TH1F *hjetpt = (TH1F*)fHistos->FindObject("histo1");
  TH1F *hjeteta = (TH1F*)fHistos->FindObject("histo2");
  TH1F *hjetphi = (TH1F*)fHistos->FindObject("histo3");
  TH1F *hjetnjet = (TH1F*)fHistos->FindObject("histo4");

  for(  Int_t jet = 0 ; jet < fNJets ; jet++ )
    {
      if (hjetpt)   hjetpt   ->Fill ( fVectJet[jet]->pt   ) ;
      if (hjeteta)  hjeteta  ->Fill ( fVectJet[jet]->eta  ) ;
      if (hjetphi)  hjetphi  ->Fill ( fVectJet[jet]->phi  ) ;
      if (hjetnjet) hjetnjet ->Fill ( fVectJet[jet]->njet ) ;
    }

  TH1F *hjets = (TH1F*)fHistos->FindObject("histo5");
  if (hjets) hjets->Fill(fNJets);

  TH1F *hleadpart = (TH1F*)fHistos->FindObject("histo6");
  if (hleadpart) hleadpart->Fill(partleadjet);

  TProfile * hprof = (TProfile*)fHistos->FindObject("histo7");
  if (hprof) hprof->Fill(ptleadjet,partleadjet);

  TH1F *hMD = (TH1F*)fHistos->FindObject("histo8");
  for(  Int_t k = 0  ; k < partleadjet ; k++)
    { hMD->Fill( zpartljet[k] ); }

  TProfile * hphi = (TProfile*)fHistos->FindObject("histo9");
  for(  Int_t k = 0  ; k < partleadjet ; k++)
    { hphi->Fill( TMath::RadToDeg() * dphipartljet [k] , fNPart ) ; }

  TProfile * htpd = (TProfile*)fHistos->FindObject("histo10");
  for(  Int_t k = 0  ; k < fNPart ; k++)
    { htpd->Fill( TMath::RadToDeg() * dphipartljet [k] , ptsumevent ) ; }


  TProfile * hprof1 = (TProfile*)fHistos->FindObject("histo21");
  if (hprof1) hprof1->Fill(ptleadjet, fNPart);

  TProfile * hprof2 = (TProfile*)fHistos->FindObject("histo22");
  if (hprof2) hprof2->Fill(ptleadjet, ptsumevent);

  TProfile * hprof1toward = (TProfile*)fHistos->FindObject("histo21_toward");
  TProfile * hprof1transverse = (TProfile*)fHistos->FindObject("histo21_transverse");
  TProfile * hprof1away = (TProfile*)fHistos->FindObject("histo21_away");
  TProfile * hprof2toward = (TProfile*)fHistos->FindObject("histo22_toward");
  TProfile * hprof2transverse = (TProfile*)fHistos->FindObject("histo22_transverse");
  TProfile * hprof2away = (TProfile*)fHistos->FindObject("histo22_away");
  TH1F * hpttoward = (TH1F*)fHistos->FindObject("histo23_toward");
  TH1F * hpttransverse = (TH1F*)fHistos->FindObject("histo23_transverse");
  TH1F * hptaway = (TH1F*)fHistos->FindObject("histo23_away");

  if ( (hprof1toward) && (hprof1transverse) && (hprof1away) && (hprof2toward) && (hprof2transverse) && (hprof2away) )
    {
      for( Int_t part = 0  ; part < fNPart ; part++)
	{
	  Double_t ptpart = fVectParticle[part]->pt ; // pt of particle
	  if ( ( dphipartljet[part] >=0.) && ( dphipartljet[part] < kPI/3. ) )
	    {
	      hprof1toward->Fill( ptleadjet, fNPart );
	      hprof2toward->Fill( ptleadjet, ptsumevent);
	      hpttoward->Fill( ptpart );
	    }
	  else
	    if ( ( dphipartljet[part] >= (kPI/3.)) && ( dphipartljet[part] < (2.*kPI/3.)) )
	      {
		hprof1transverse->Fill( ptleadjet, fNPart );
		hprof2transverse->Fill( ptleadjet, ptsumevent);
		hpttransverse->Fill( ptpart );
	      }
	    else
	      if ( ( dphipartljet[part] >= ( 2.*kPI/3.)) && ( dphipartljet[part] < kPI ) )
		{
		  hprof1away->Fill( ptleadjet, fNPart );
		  hprof2away->Fill( ptleadjet, ptsumevent);
		  hptaway->Fill( ptpart );
		}
	}
    }

  delete [] dphipartljet;
  delete [] zpartljet;
  delete [] idxpartLJ;

}

//______________________________________________________________________________
void AliCdfJetFinder::Clean()
{
  // CLEANING SECTION
  for (  Int_t i = 0 ; i < fNPart ; i++ ){
    delete fVectParticle[i];
    fVectParticle[i] = 0;
  }
  delete [] fVectParticle;fVectParticle = 0;

  for (  Int_t i = 0 ; i < fNJets ; i++ ){
    delete fVectJet[i];
    fVectJet[i] = 0;
  } 
  delete [] fVectJet;fVectJet = 0;

  delete [] fPtArray;fPtArray = 0;
  delete [] fIdxArray;fIdxArray = 0;

  Reset();

}

 AliCdfJetFinder.cxx:1
 AliCdfJetFinder.cxx:2
 AliCdfJetFinder.cxx:3
 AliCdfJetFinder.cxx:4
 AliCdfJetFinder.cxx:5
 AliCdfJetFinder.cxx:6
 AliCdfJetFinder.cxx:7
 AliCdfJetFinder.cxx:8
 AliCdfJetFinder.cxx:9
 AliCdfJetFinder.cxx:10
 AliCdfJetFinder.cxx:11
 AliCdfJetFinder.cxx:12
 AliCdfJetFinder.cxx:13
 AliCdfJetFinder.cxx:14
 AliCdfJetFinder.cxx:15
 AliCdfJetFinder.cxx:16
 AliCdfJetFinder.cxx:17
 AliCdfJetFinder.cxx:18
 AliCdfJetFinder.cxx:19
 AliCdfJetFinder.cxx:20
 AliCdfJetFinder.cxx:21
 AliCdfJetFinder.cxx:22
 AliCdfJetFinder.cxx:23
 AliCdfJetFinder.cxx:24
 AliCdfJetFinder.cxx:25
 AliCdfJetFinder.cxx:26
 AliCdfJetFinder.cxx:27
 AliCdfJetFinder.cxx:28
 AliCdfJetFinder.cxx:29
 AliCdfJetFinder.cxx:30
 AliCdfJetFinder.cxx:31
 AliCdfJetFinder.cxx:32
 AliCdfJetFinder.cxx:33
 AliCdfJetFinder.cxx:34
 AliCdfJetFinder.cxx:35
 AliCdfJetFinder.cxx:36
 AliCdfJetFinder.cxx:37
 AliCdfJetFinder.cxx:38
 AliCdfJetFinder.cxx:39
 AliCdfJetFinder.cxx:40
 AliCdfJetFinder.cxx:41
 AliCdfJetFinder.cxx:42
 AliCdfJetFinder.cxx:43
 AliCdfJetFinder.cxx:44
 AliCdfJetFinder.cxx:45
 AliCdfJetFinder.cxx:46
 AliCdfJetFinder.cxx:47
 AliCdfJetFinder.cxx:48
 AliCdfJetFinder.cxx:49
 AliCdfJetFinder.cxx:50
 AliCdfJetFinder.cxx:51
 AliCdfJetFinder.cxx:52
 AliCdfJetFinder.cxx:53
 AliCdfJetFinder.cxx:54
 AliCdfJetFinder.cxx:55
 AliCdfJetFinder.cxx:56
 AliCdfJetFinder.cxx:57
 AliCdfJetFinder.cxx:58
 AliCdfJetFinder.cxx:59
 AliCdfJetFinder.cxx:60
 AliCdfJetFinder.cxx:61
 AliCdfJetFinder.cxx:62
 AliCdfJetFinder.cxx:63
 AliCdfJetFinder.cxx:64
 AliCdfJetFinder.cxx:65
 AliCdfJetFinder.cxx:66
 AliCdfJetFinder.cxx:67
 AliCdfJetFinder.cxx:68
 AliCdfJetFinder.cxx:69
 AliCdfJetFinder.cxx:70
 AliCdfJetFinder.cxx:71
 AliCdfJetFinder.cxx:72
 AliCdfJetFinder.cxx:73
 AliCdfJetFinder.cxx:74
 AliCdfJetFinder.cxx:75
 AliCdfJetFinder.cxx:76
 AliCdfJetFinder.cxx:77
 AliCdfJetFinder.cxx:78
 AliCdfJetFinder.cxx:79
 AliCdfJetFinder.cxx:80
 AliCdfJetFinder.cxx:81
 AliCdfJetFinder.cxx:82
 AliCdfJetFinder.cxx:83
 AliCdfJetFinder.cxx:84
 AliCdfJetFinder.cxx:85
 AliCdfJetFinder.cxx:86
 AliCdfJetFinder.cxx:87
 AliCdfJetFinder.cxx:88
 AliCdfJetFinder.cxx:89
 AliCdfJetFinder.cxx:90
 AliCdfJetFinder.cxx:91
 AliCdfJetFinder.cxx:92
 AliCdfJetFinder.cxx:93
 AliCdfJetFinder.cxx:94
 AliCdfJetFinder.cxx:95
 AliCdfJetFinder.cxx:96
 AliCdfJetFinder.cxx:97
 AliCdfJetFinder.cxx:98
 AliCdfJetFinder.cxx:99
 AliCdfJetFinder.cxx:100
 AliCdfJetFinder.cxx:101
 AliCdfJetFinder.cxx:102
 AliCdfJetFinder.cxx:103
 AliCdfJetFinder.cxx:104
 AliCdfJetFinder.cxx:105
 AliCdfJetFinder.cxx:106
 AliCdfJetFinder.cxx:107
 AliCdfJetFinder.cxx:108
 AliCdfJetFinder.cxx:109
 AliCdfJetFinder.cxx:110
 AliCdfJetFinder.cxx:111
 AliCdfJetFinder.cxx:112
 AliCdfJetFinder.cxx:113
 AliCdfJetFinder.cxx:114
 AliCdfJetFinder.cxx:115
 AliCdfJetFinder.cxx:116
 AliCdfJetFinder.cxx:117
 AliCdfJetFinder.cxx:118
 AliCdfJetFinder.cxx:119
 AliCdfJetFinder.cxx:120
 AliCdfJetFinder.cxx:121
 AliCdfJetFinder.cxx:122
 AliCdfJetFinder.cxx:123
 AliCdfJetFinder.cxx:124
 AliCdfJetFinder.cxx:125
 AliCdfJetFinder.cxx:126
 AliCdfJetFinder.cxx:127
 AliCdfJetFinder.cxx:128
 AliCdfJetFinder.cxx:129
 AliCdfJetFinder.cxx:130
 AliCdfJetFinder.cxx:131
 AliCdfJetFinder.cxx:132
 AliCdfJetFinder.cxx:133
 AliCdfJetFinder.cxx:134
 AliCdfJetFinder.cxx:135
 AliCdfJetFinder.cxx:136
 AliCdfJetFinder.cxx:137
 AliCdfJetFinder.cxx:138
 AliCdfJetFinder.cxx:139
 AliCdfJetFinder.cxx:140
 AliCdfJetFinder.cxx:141
 AliCdfJetFinder.cxx:142
 AliCdfJetFinder.cxx:143
 AliCdfJetFinder.cxx:144
 AliCdfJetFinder.cxx:145
 AliCdfJetFinder.cxx:146
 AliCdfJetFinder.cxx:147
 AliCdfJetFinder.cxx:148
 AliCdfJetFinder.cxx:149
 AliCdfJetFinder.cxx:150
 AliCdfJetFinder.cxx:151
 AliCdfJetFinder.cxx:152
 AliCdfJetFinder.cxx:153
 AliCdfJetFinder.cxx:154
 AliCdfJetFinder.cxx:155
 AliCdfJetFinder.cxx:156
 AliCdfJetFinder.cxx:157
 AliCdfJetFinder.cxx:158
 AliCdfJetFinder.cxx:159
 AliCdfJetFinder.cxx:160
 AliCdfJetFinder.cxx:161
 AliCdfJetFinder.cxx:162
 AliCdfJetFinder.cxx:163
 AliCdfJetFinder.cxx:164
 AliCdfJetFinder.cxx:165
 AliCdfJetFinder.cxx:166
 AliCdfJetFinder.cxx:167
 AliCdfJetFinder.cxx:168
 AliCdfJetFinder.cxx:169
 AliCdfJetFinder.cxx:170
 AliCdfJetFinder.cxx:171
 AliCdfJetFinder.cxx:172
 AliCdfJetFinder.cxx:173
 AliCdfJetFinder.cxx:174
 AliCdfJetFinder.cxx:175
 AliCdfJetFinder.cxx:176
 AliCdfJetFinder.cxx:177
 AliCdfJetFinder.cxx:178
 AliCdfJetFinder.cxx:179
 AliCdfJetFinder.cxx:180
 AliCdfJetFinder.cxx:181
 AliCdfJetFinder.cxx:182
 AliCdfJetFinder.cxx:183
 AliCdfJetFinder.cxx:184
 AliCdfJetFinder.cxx:185
 AliCdfJetFinder.cxx:186
 AliCdfJetFinder.cxx:187
 AliCdfJetFinder.cxx:188
 AliCdfJetFinder.cxx:189
 AliCdfJetFinder.cxx:190
 AliCdfJetFinder.cxx:191
 AliCdfJetFinder.cxx:192
 AliCdfJetFinder.cxx:193
 AliCdfJetFinder.cxx:194
 AliCdfJetFinder.cxx:195
 AliCdfJetFinder.cxx:196
 AliCdfJetFinder.cxx:197
 AliCdfJetFinder.cxx:198
 AliCdfJetFinder.cxx:199
 AliCdfJetFinder.cxx:200
 AliCdfJetFinder.cxx:201
 AliCdfJetFinder.cxx:202
 AliCdfJetFinder.cxx:203
 AliCdfJetFinder.cxx:204
 AliCdfJetFinder.cxx:205
 AliCdfJetFinder.cxx:206
 AliCdfJetFinder.cxx:207
 AliCdfJetFinder.cxx:208
 AliCdfJetFinder.cxx:209
 AliCdfJetFinder.cxx:210
 AliCdfJetFinder.cxx:211
 AliCdfJetFinder.cxx:212
 AliCdfJetFinder.cxx:213
 AliCdfJetFinder.cxx:214
 AliCdfJetFinder.cxx:215
 AliCdfJetFinder.cxx:216
 AliCdfJetFinder.cxx:217
 AliCdfJetFinder.cxx:218
 AliCdfJetFinder.cxx:219
 AliCdfJetFinder.cxx:220
 AliCdfJetFinder.cxx:221
 AliCdfJetFinder.cxx:222
 AliCdfJetFinder.cxx:223
 AliCdfJetFinder.cxx:224
 AliCdfJetFinder.cxx:225
 AliCdfJetFinder.cxx:226
 AliCdfJetFinder.cxx:227
 AliCdfJetFinder.cxx:228
 AliCdfJetFinder.cxx:229
 AliCdfJetFinder.cxx:230
 AliCdfJetFinder.cxx:231
 AliCdfJetFinder.cxx:232
 AliCdfJetFinder.cxx:233
 AliCdfJetFinder.cxx:234
 AliCdfJetFinder.cxx:235
 AliCdfJetFinder.cxx:236
 AliCdfJetFinder.cxx:237
 AliCdfJetFinder.cxx:238
 AliCdfJetFinder.cxx:239
 AliCdfJetFinder.cxx:240
 AliCdfJetFinder.cxx:241
 AliCdfJetFinder.cxx:242
 AliCdfJetFinder.cxx:243
 AliCdfJetFinder.cxx:244
 AliCdfJetFinder.cxx:245
 AliCdfJetFinder.cxx:246
 AliCdfJetFinder.cxx:247
 AliCdfJetFinder.cxx:248
 AliCdfJetFinder.cxx:249
 AliCdfJetFinder.cxx:250
 AliCdfJetFinder.cxx:251
 AliCdfJetFinder.cxx:252
 AliCdfJetFinder.cxx:253
 AliCdfJetFinder.cxx:254
 AliCdfJetFinder.cxx:255
 AliCdfJetFinder.cxx:256
 AliCdfJetFinder.cxx:257
 AliCdfJetFinder.cxx:258
 AliCdfJetFinder.cxx:259
 AliCdfJetFinder.cxx:260
 AliCdfJetFinder.cxx:261
 AliCdfJetFinder.cxx:262
 AliCdfJetFinder.cxx:263
 AliCdfJetFinder.cxx:264
 AliCdfJetFinder.cxx:265
 AliCdfJetFinder.cxx:266
 AliCdfJetFinder.cxx:267
 AliCdfJetFinder.cxx:268
 AliCdfJetFinder.cxx:269
 AliCdfJetFinder.cxx:270
 AliCdfJetFinder.cxx:271
 AliCdfJetFinder.cxx:272
 AliCdfJetFinder.cxx:273
 AliCdfJetFinder.cxx:274
 AliCdfJetFinder.cxx:275
 AliCdfJetFinder.cxx:276
 AliCdfJetFinder.cxx:277
 AliCdfJetFinder.cxx:278
 AliCdfJetFinder.cxx:279
 AliCdfJetFinder.cxx:280
 AliCdfJetFinder.cxx:281
 AliCdfJetFinder.cxx:282
 AliCdfJetFinder.cxx:283
 AliCdfJetFinder.cxx:284
 AliCdfJetFinder.cxx:285
 AliCdfJetFinder.cxx:286
 AliCdfJetFinder.cxx:287
 AliCdfJetFinder.cxx:288
 AliCdfJetFinder.cxx:289
 AliCdfJetFinder.cxx:290
 AliCdfJetFinder.cxx:291
 AliCdfJetFinder.cxx:292
 AliCdfJetFinder.cxx:293
 AliCdfJetFinder.cxx:294
 AliCdfJetFinder.cxx:295
 AliCdfJetFinder.cxx:296
 AliCdfJetFinder.cxx:297
 AliCdfJetFinder.cxx:298
 AliCdfJetFinder.cxx:299
 AliCdfJetFinder.cxx:300
 AliCdfJetFinder.cxx:301
 AliCdfJetFinder.cxx:302
 AliCdfJetFinder.cxx:303
 AliCdfJetFinder.cxx:304
 AliCdfJetFinder.cxx:305
 AliCdfJetFinder.cxx:306
 AliCdfJetFinder.cxx:307
 AliCdfJetFinder.cxx:308
 AliCdfJetFinder.cxx:309
 AliCdfJetFinder.cxx:310
 AliCdfJetFinder.cxx:311
 AliCdfJetFinder.cxx:312
 AliCdfJetFinder.cxx:313
 AliCdfJetFinder.cxx:314
 AliCdfJetFinder.cxx:315
 AliCdfJetFinder.cxx:316
 AliCdfJetFinder.cxx:317
 AliCdfJetFinder.cxx:318
 AliCdfJetFinder.cxx:319
 AliCdfJetFinder.cxx:320
 AliCdfJetFinder.cxx:321
 AliCdfJetFinder.cxx:322
 AliCdfJetFinder.cxx:323
 AliCdfJetFinder.cxx:324
 AliCdfJetFinder.cxx:325
 AliCdfJetFinder.cxx:326
 AliCdfJetFinder.cxx:327
 AliCdfJetFinder.cxx:328
 AliCdfJetFinder.cxx:329
 AliCdfJetFinder.cxx:330
 AliCdfJetFinder.cxx:331
 AliCdfJetFinder.cxx:332
 AliCdfJetFinder.cxx:333
 AliCdfJetFinder.cxx:334
 AliCdfJetFinder.cxx:335
 AliCdfJetFinder.cxx:336
 AliCdfJetFinder.cxx:337
 AliCdfJetFinder.cxx:338
 AliCdfJetFinder.cxx:339
 AliCdfJetFinder.cxx:340
 AliCdfJetFinder.cxx:341
 AliCdfJetFinder.cxx:342
 AliCdfJetFinder.cxx:343
 AliCdfJetFinder.cxx:344
 AliCdfJetFinder.cxx:345
 AliCdfJetFinder.cxx:346
 AliCdfJetFinder.cxx:347
 AliCdfJetFinder.cxx:348
 AliCdfJetFinder.cxx:349
 AliCdfJetFinder.cxx:350
 AliCdfJetFinder.cxx:351
 AliCdfJetFinder.cxx:352
 AliCdfJetFinder.cxx:353
 AliCdfJetFinder.cxx:354
 AliCdfJetFinder.cxx:355
 AliCdfJetFinder.cxx:356
 AliCdfJetFinder.cxx:357
 AliCdfJetFinder.cxx:358
 AliCdfJetFinder.cxx:359
 AliCdfJetFinder.cxx:360
 AliCdfJetFinder.cxx:361
 AliCdfJetFinder.cxx:362
 AliCdfJetFinder.cxx:363
 AliCdfJetFinder.cxx:364
 AliCdfJetFinder.cxx:365
 AliCdfJetFinder.cxx:366
 AliCdfJetFinder.cxx:367
 AliCdfJetFinder.cxx:368
 AliCdfJetFinder.cxx:369
 AliCdfJetFinder.cxx:370
 AliCdfJetFinder.cxx:371
 AliCdfJetFinder.cxx:372
 AliCdfJetFinder.cxx:373
 AliCdfJetFinder.cxx:374
 AliCdfJetFinder.cxx:375
 AliCdfJetFinder.cxx:376
 AliCdfJetFinder.cxx:377
 AliCdfJetFinder.cxx:378
 AliCdfJetFinder.cxx:379
 AliCdfJetFinder.cxx:380
 AliCdfJetFinder.cxx:381
 AliCdfJetFinder.cxx:382
 AliCdfJetFinder.cxx:383
 AliCdfJetFinder.cxx:384
 AliCdfJetFinder.cxx:385
 AliCdfJetFinder.cxx:386
 AliCdfJetFinder.cxx:387
 AliCdfJetFinder.cxx:388
 AliCdfJetFinder.cxx:389
 AliCdfJetFinder.cxx:390
 AliCdfJetFinder.cxx:391
 AliCdfJetFinder.cxx:392
 AliCdfJetFinder.cxx:393
 AliCdfJetFinder.cxx:394
 AliCdfJetFinder.cxx:395
 AliCdfJetFinder.cxx:396
 AliCdfJetFinder.cxx:397
 AliCdfJetFinder.cxx:398
 AliCdfJetFinder.cxx:399
 AliCdfJetFinder.cxx:400
 AliCdfJetFinder.cxx:401
 AliCdfJetFinder.cxx:402
 AliCdfJetFinder.cxx:403
 AliCdfJetFinder.cxx:404
 AliCdfJetFinder.cxx:405
 AliCdfJetFinder.cxx:406
 AliCdfJetFinder.cxx:407
 AliCdfJetFinder.cxx:408
 AliCdfJetFinder.cxx:409
 AliCdfJetFinder.cxx:410
 AliCdfJetFinder.cxx:411
 AliCdfJetFinder.cxx:412
 AliCdfJetFinder.cxx:413
 AliCdfJetFinder.cxx:414
 AliCdfJetFinder.cxx:415
 AliCdfJetFinder.cxx:416
 AliCdfJetFinder.cxx:417
 AliCdfJetFinder.cxx:418
 AliCdfJetFinder.cxx:419
 AliCdfJetFinder.cxx:420
 AliCdfJetFinder.cxx:421
 AliCdfJetFinder.cxx:422
 AliCdfJetFinder.cxx:423
 AliCdfJetFinder.cxx:424
 AliCdfJetFinder.cxx:425
 AliCdfJetFinder.cxx:426
 AliCdfJetFinder.cxx:427
 AliCdfJetFinder.cxx:428
 AliCdfJetFinder.cxx:429
 AliCdfJetFinder.cxx:430
 AliCdfJetFinder.cxx:431
 AliCdfJetFinder.cxx:432
 AliCdfJetFinder.cxx:433
 AliCdfJetFinder.cxx:434
 AliCdfJetFinder.cxx:435
 AliCdfJetFinder.cxx:436
 AliCdfJetFinder.cxx:437
 AliCdfJetFinder.cxx:438
 AliCdfJetFinder.cxx:439
 AliCdfJetFinder.cxx:440
 AliCdfJetFinder.cxx:441
 AliCdfJetFinder.cxx:442
 AliCdfJetFinder.cxx:443
 AliCdfJetFinder.cxx:444
 AliCdfJetFinder.cxx:445
 AliCdfJetFinder.cxx:446
 AliCdfJetFinder.cxx:447
 AliCdfJetFinder.cxx:448
 AliCdfJetFinder.cxx:449
 AliCdfJetFinder.cxx:450
 AliCdfJetFinder.cxx:451
 AliCdfJetFinder.cxx:452
 AliCdfJetFinder.cxx:453
 AliCdfJetFinder.cxx:454
 AliCdfJetFinder.cxx:455
 AliCdfJetFinder.cxx:456
 AliCdfJetFinder.cxx:457
 AliCdfJetFinder.cxx:458
 AliCdfJetFinder.cxx:459
 AliCdfJetFinder.cxx:460
 AliCdfJetFinder.cxx:461
 AliCdfJetFinder.cxx:462
 AliCdfJetFinder.cxx:463
 AliCdfJetFinder.cxx:464
 AliCdfJetFinder.cxx:465
 AliCdfJetFinder.cxx:466
 AliCdfJetFinder.cxx:467
 AliCdfJetFinder.cxx:468
 AliCdfJetFinder.cxx:469
 AliCdfJetFinder.cxx:470
 AliCdfJetFinder.cxx:471
 AliCdfJetFinder.cxx:472
 AliCdfJetFinder.cxx:473
 AliCdfJetFinder.cxx:474
 AliCdfJetFinder.cxx:475
 AliCdfJetFinder.cxx:476
 AliCdfJetFinder.cxx:477
 AliCdfJetFinder.cxx:478
 AliCdfJetFinder.cxx:479
 AliCdfJetFinder.cxx:480
 AliCdfJetFinder.cxx:481
 AliCdfJetFinder.cxx:482
 AliCdfJetFinder.cxx:483
 AliCdfJetFinder.cxx:484
 AliCdfJetFinder.cxx:485
 AliCdfJetFinder.cxx:486
 AliCdfJetFinder.cxx:487
 AliCdfJetFinder.cxx:488
 AliCdfJetFinder.cxx:489
 AliCdfJetFinder.cxx:490
 AliCdfJetFinder.cxx:491
 AliCdfJetFinder.cxx:492
 AliCdfJetFinder.cxx:493
 AliCdfJetFinder.cxx:494
 AliCdfJetFinder.cxx:495
 AliCdfJetFinder.cxx:496
 AliCdfJetFinder.cxx:497
 AliCdfJetFinder.cxx:498
 AliCdfJetFinder.cxx:499
 AliCdfJetFinder.cxx:500
 AliCdfJetFinder.cxx:501
 AliCdfJetFinder.cxx:502
 AliCdfJetFinder.cxx:503
 AliCdfJetFinder.cxx:504
 AliCdfJetFinder.cxx:505
 AliCdfJetFinder.cxx:506
 AliCdfJetFinder.cxx:507
 AliCdfJetFinder.cxx:508
 AliCdfJetFinder.cxx:509
 AliCdfJetFinder.cxx:510
 AliCdfJetFinder.cxx:511
 AliCdfJetFinder.cxx:512
 AliCdfJetFinder.cxx:513
 AliCdfJetFinder.cxx:514
 AliCdfJetFinder.cxx:515
 AliCdfJetFinder.cxx:516
 AliCdfJetFinder.cxx:517
 AliCdfJetFinder.cxx:518
 AliCdfJetFinder.cxx:519
 AliCdfJetFinder.cxx:520
 AliCdfJetFinder.cxx:521
 AliCdfJetFinder.cxx:522
 AliCdfJetFinder.cxx:523
 AliCdfJetFinder.cxx:524
 AliCdfJetFinder.cxx:525
 AliCdfJetFinder.cxx:526
 AliCdfJetFinder.cxx:527
 AliCdfJetFinder.cxx:528
 AliCdfJetFinder.cxx:529
 AliCdfJetFinder.cxx:530
 AliCdfJetFinder.cxx:531
 AliCdfJetFinder.cxx:532
 AliCdfJetFinder.cxx:533
 AliCdfJetFinder.cxx:534
 AliCdfJetFinder.cxx:535
 AliCdfJetFinder.cxx:536
 AliCdfJetFinder.cxx:537
 AliCdfJetFinder.cxx:538
 AliCdfJetFinder.cxx:539
 AliCdfJetFinder.cxx:540
 AliCdfJetFinder.cxx:541
 AliCdfJetFinder.cxx:542
 AliCdfJetFinder.cxx:543
 AliCdfJetFinder.cxx:544
 AliCdfJetFinder.cxx:545
 AliCdfJetFinder.cxx:546
 AliCdfJetFinder.cxx:547
 AliCdfJetFinder.cxx:548
 AliCdfJetFinder.cxx:549
 AliCdfJetFinder.cxx:550
 AliCdfJetFinder.cxx:551
 AliCdfJetFinder.cxx:552
 AliCdfJetFinder.cxx:553
 AliCdfJetFinder.cxx:554
 AliCdfJetFinder.cxx:555
 AliCdfJetFinder.cxx:556
 AliCdfJetFinder.cxx:557
 AliCdfJetFinder.cxx:558
 AliCdfJetFinder.cxx:559
 AliCdfJetFinder.cxx:560
 AliCdfJetFinder.cxx:561
 AliCdfJetFinder.cxx:562
 AliCdfJetFinder.cxx:563
 AliCdfJetFinder.cxx:564
 AliCdfJetFinder.cxx:565
 AliCdfJetFinder.cxx:566
 AliCdfJetFinder.cxx:567
 AliCdfJetFinder.cxx:568
 AliCdfJetFinder.cxx:569
 AliCdfJetFinder.cxx:570
 AliCdfJetFinder.cxx:571
 AliCdfJetFinder.cxx:572
 AliCdfJetFinder.cxx:573
 AliCdfJetFinder.cxx:574
 AliCdfJetFinder.cxx:575
 AliCdfJetFinder.cxx:576
 AliCdfJetFinder.cxx:577
 AliCdfJetFinder.cxx:578
 AliCdfJetFinder.cxx:579
 AliCdfJetFinder.cxx:580
 AliCdfJetFinder.cxx:581
 AliCdfJetFinder.cxx:582
 AliCdfJetFinder.cxx:583
 AliCdfJetFinder.cxx:584
 AliCdfJetFinder.cxx:585
 AliCdfJetFinder.cxx:586
 AliCdfJetFinder.cxx:587
 AliCdfJetFinder.cxx:588
 AliCdfJetFinder.cxx:589
 AliCdfJetFinder.cxx:590
 AliCdfJetFinder.cxx:591
 AliCdfJetFinder.cxx:592
 AliCdfJetFinder.cxx:593
 AliCdfJetFinder.cxx:594
 AliCdfJetFinder.cxx:595
 AliCdfJetFinder.cxx:596
 AliCdfJetFinder.cxx:597
 AliCdfJetFinder.cxx:598
 AliCdfJetFinder.cxx:599
 AliCdfJetFinder.cxx:600
 AliCdfJetFinder.cxx:601
 AliCdfJetFinder.cxx:602
 AliCdfJetFinder.cxx:603
 AliCdfJetFinder.cxx:604
 AliCdfJetFinder.cxx:605
 AliCdfJetFinder.cxx:606
 AliCdfJetFinder.cxx:607
 AliCdfJetFinder.cxx:608
 AliCdfJetFinder.cxx:609
 AliCdfJetFinder.cxx:610
 AliCdfJetFinder.cxx:611
 AliCdfJetFinder.cxx:612
 AliCdfJetFinder.cxx:613
 AliCdfJetFinder.cxx:614
 AliCdfJetFinder.cxx:615
 AliCdfJetFinder.cxx:616
 AliCdfJetFinder.cxx:617
 AliCdfJetFinder.cxx:618
 AliCdfJetFinder.cxx:619
 AliCdfJetFinder.cxx:620
 AliCdfJetFinder.cxx:621
 AliCdfJetFinder.cxx:622
 AliCdfJetFinder.cxx:623
 AliCdfJetFinder.cxx:624
 AliCdfJetFinder.cxx:625
 AliCdfJetFinder.cxx:626
 AliCdfJetFinder.cxx:627
 AliCdfJetFinder.cxx:628
 AliCdfJetFinder.cxx:629
 AliCdfJetFinder.cxx:630
 AliCdfJetFinder.cxx:631
 AliCdfJetFinder.cxx:632
 AliCdfJetFinder.cxx:633
 AliCdfJetFinder.cxx:634
 AliCdfJetFinder.cxx:635
 AliCdfJetFinder.cxx:636
 AliCdfJetFinder.cxx:637
 AliCdfJetFinder.cxx:638
 AliCdfJetFinder.cxx:639
 AliCdfJetFinder.cxx:640
 AliCdfJetFinder.cxx:641
 AliCdfJetFinder.cxx:642
 AliCdfJetFinder.cxx:643
 AliCdfJetFinder.cxx:644
 AliCdfJetFinder.cxx:645
 AliCdfJetFinder.cxx:646
 AliCdfJetFinder.cxx:647
 AliCdfJetFinder.cxx:648
 AliCdfJetFinder.cxx:649
 AliCdfJetFinder.cxx:650
 AliCdfJetFinder.cxx:651
 AliCdfJetFinder.cxx:652
 AliCdfJetFinder.cxx:653
 AliCdfJetFinder.cxx:654
 AliCdfJetFinder.cxx:655
 AliCdfJetFinder.cxx:656
 AliCdfJetFinder.cxx:657
 AliCdfJetFinder.cxx:658
 AliCdfJetFinder.cxx:659
 AliCdfJetFinder.cxx:660
 AliCdfJetFinder.cxx:661
 AliCdfJetFinder.cxx:662
 AliCdfJetFinder.cxx:663
 AliCdfJetFinder.cxx:664
 AliCdfJetFinder.cxx:665
 AliCdfJetFinder.cxx:666
 AliCdfJetFinder.cxx:667
 AliCdfJetFinder.cxx:668
 AliCdfJetFinder.cxx:669
 AliCdfJetFinder.cxx:670
 AliCdfJetFinder.cxx:671
 AliCdfJetFinder.cxx:672
 AliCdfJetFinder.cxx:673
 AliCdfJetFinder.cxx:674
 AliCdfJetFinder.cxx:675
 AliCdfJetFinder.cxx:676
 AliCdfJetFinder.cxx:677
 AliCdfJetFinder.cxx:678
 AliCdfJetFinder.cxx:679
 AliCdfJetFinder.cxx:680
 AliCdfJetFinder.cxx:681
 AliCdfJetFinder.cxx:682
 AliCdfJetFinder.cxx:683
 AliCdfJetFinder.cxx:684
 AliCdfJetFinder.cxx:685
 AliCdfJetFinder.cxx:686
 AliCdfJetFinder.cxx:687
 AliCdfJetFinder.cxx:688
 AliCdfJetFinder.cxx:689
 AliCdfJetFinder.cxx:690
 AliCdfJetFinder.cxx:691
 AliCdfJetFinder.cxx:692
 AliCdfJetFinder.cxx:693
 AliCdfJetFinder.cxx:694
 AliCdfJetFinder.cxx:695
 AliCdfJetFinder.cxx:696
 AliCdfJetFinder.cxx:697
 AliCdfJetFinder.cxx:698
 AliCdfJetFinder.cxx:699
 AliCdfJetFinder.cxx:700
 AliCdfJetFinder.cxx:701
 AliCdfJetFinder.cxx:702
 AliCdfJetFinder.cxx:703
 AliCdfJetFinder.cxx:704
 AliCdfJetFinder.cxx:705
 AliCdfJetFinder.cxx:706
 AliCdfJetFinder.cxx:707
 AliCdfJetFinder.cxx:708
 AliCdfJetFinder.cxx:709
 AliCdfJetFinder.cxx:710
 AliCdfJetFinder.cxx:711
 AliCdfJetFinder.cxx:712
 AliCdfJetFinder.cxx:713
 AliCdfJetFinder.cxx:714
 AliCdfJetFinder.cxx:715
 AliCdfJetFinder.cxx:716
 AliCdfJetFinder.cxx:717
 AliCdfJetFinder.cxx:718
 AliCdfJetFinder.cxx:719
 AliCdfJetFinder.cxx:720
 AliCdfJetFinder.cxx:721
 AliCdfJetFinder.cxx:722
 AliCdfJetFinder.cxx:723
 AliCdfJetFinder.cxx:724
 AliCdfJetFinder.cxx:725
 AliCdfJetFinder.cxx:726
 AliCdfJetFinder.cxx:727
 AliCdfJetFinder.cxx:728
 AliCdfJetFinder.cxx:729
 AliCdfJetFinder.cxx:730
 AliCdfJetFinder.cxx:731
 AliCdfJetFinder.cxx:732
 AliCdfJetFinder.cxx:733
 AliCdfJetFinder.cxx:734
 AliCdfJetFinder.cxx:735
 AliCdfJetFinder.cxx:736
 AliCdfJetFinder.cxx:737
 AliCdfJetFinder.cxx:738
 AliCdfJetFinder.cxx:739
 AliCdfJetFinder.cxx:740
 AliCdfJetFinder.cxx:741
 AliCdfJetFinder.cxx:742
 AliCdfJetFinder.cxx:743
 AliCdfJetFinder.cxx:744
 AliCdfJetFinder.cxx:745
 AliCdfJetFinder.cxx:746
 AliCdfJetFinder.cxx:747
 AliCdfJetFinder.cxx:748
 AliCdfJetFinder.cxx:749
 AliCdfJetFinder.cxx:750
 AliCdfJetFinder.cxx:751
 AliCdfJetFinder.cxx:752
 AliCdfJetFinder.cxx:753
 AliCdfJetFinder.cxx:754
 AliCdfJetFinder.cxx:755
 AliCdfJetFinder.cxx:756
 AliCdfJetFinder.cxx:757
 AliCdfJetFinder.cxx:758
 AliCdfJetFinder.cxx:759
 AliCdfJetFinder.cxx:760
 AliCdfJetFinder.cxx:761
 AliCdfJetFinder.cxx:762
 AliCdfJetFinder.cxx:763
 AliCdfJetFinder.cxx:764
 AliCdfJetFinder.cxx:765
 AliCdfJetFinder.cxx:766
 AliCdfJetFinder.cxx:767
 AliCdfJetFinder.cxx:768
 AliCdfJetFinder.cxx:769
 AliCdfJetFinder.cxx:770
 AliCdfJetFinder.cxx:771
 AliCdfJetFinder.cxx:772
 AliCdfJetFinder.cxx:773
 AliCdfJetFinder.cxx:774
 AliCdfJetFinder.cxx:775
 AliCdfJetFinder.cxx:776
 AliCdfJetFinder.cxx:777
 AliCdfJetFinder.cxx:778
 AliCdfJetFinder.cxx:779
 AliCdfJetFinder.cxx:780
 AliCdfJetFinder.cxx:781
 AliCdfJetFinder.cxx:782
 AliCdfJetFinder.cxx:783
 AliCdfJetFinder.cxx:784
 AliCdfJetFinder.cxx:785
 AliCdfJetFinder.cxx:786
 AliCdfJetFinder.cxx:787
 AliCdfJetFinder.cxx:788
 AliCdfJetFinder.cxx:789
 AliCdfJetFinder.cxx:790
 AliCdfJetFinder.cxx:791
 AliCdfJetFinder.cxx:792
 AliCdfJetFinder.cxx:793
 AliCdfJetFinder.cxx:794
 AliCdfJetFinder.cxx:795
 AliCdfJetFinder.cxx:796
 AliCdfJetFinder.cxx:797
 AliCdfJetFinder.cxx:798
 AliCdfJetFinder.cxx:799
 AliCdfJetFinder.cxx:800
 AliCdfJetFinder.cxx:801
 AliCdfJetFinder.cxx:802
 AliCdfJetFinder.cxx:803
 AliCdfJetFinder.cxx:804
 AliCdfJetFinder.cxx:805
 AliCdfJetFinder.cxx:806
 AliCdfJetFinder.cxx:807
 AliCdfJetFinder.cxx:808
 AliCdfJetFinder.cxx:809
 AliCdfJetFinder.cxx:810
 AliCdfJetFinder.cxx:811
 AliCdfJetFinder.cxx:812
 AliCdfJetFinder.cxx:813
 AliCdfJetFinder.cxx:814
 AliCdfJetFinder.cxx:815
 AliCdfJetFinder.cxx:816
 AliCdfJetFinder.cxx:817
 AliCdfJetFinder.cxx:818
 AliCdfJetFinder.cxx:819
 AliCdfJetFinder.cxx:820
 AliCdfJetFinder.cxx:821
 AliCdfJetFinder.cxx:822
 AliCdfJetFinder.cxx:823
 AliCdfJetFinder.cxx:824
 AliCdfJetFinder.cxx:825
 AliCdfJetFinder.cxx:826
 AliCdfJetFinder.cxx:827
 AliCdfJetFinder.cxx:828
 AliCdfJetFinder.cxx:829
 AliCdfJetFinder.cxx:830
 AliCdfJetFinder.cxx:831
 AliCdfJetFinder.cxx:832
 AliCdfJetFinder.cxx:833
 AliCdfJetFinder.cxx:834
 AliCdfJetFinder.cxx:835
 AliCdfJetFinder.cxx:836
 AliCdfJetFinder.cxx:837
 AliCdfJetFinder.cxx:838
 AliCdfJetFinder.cxx:839
 AliCdfJetFinder.cxx:840
 AliCdfJetFinder.cxx:841
 AliCdfJetFinder.cxx:842
 AliCdfJetFinder.cxx:843
 AliCdfJetFinder.cxx:844
 AliCdfJetFinder.cxx:845
 AliCdfJetFinder.cxx:846
 AliCdfJetFinder.cxx:847
 AliCdfJetFinder.cxx:848
 AliCdfJetFinder.cxx:849
 AliCdfJetFinder.cxx:850
 AliCdfJetFinder.cxx:851
 AliCdfJetFinder.cxx:852
 AliCdfJetFinder.cxx:853
 AliCdfJetFinder.cxx:854
 AliCdfJetFinder.cxx:855
 AliCdfJetFinder.cxx:856
 AliCdfJetFinder.cxx:857
 AliCdfJetFinder.cxx:858
 AliCdfJetFinder.cxx:859
 AliCdfJetFinder.cxx:860
 AliCdfJetFinder.cxx:861
 AliCdfJetFinder.cxx:862
 AliCdfJetFinder.cxx:863
 AliCdfJetFinder.cxx:864
 AliCdfJetFinder.cxx:865
 AliCdfJetFinder.cxx:866
 AliCdfJetFinder.cxx:867
 AliCdfJetFinder.cxx:868
 AliCdfJetFinder.cxx:869
 AliCdfJetFinder.cxx:870
 AliCdfJetFinder.cxx:871
 AliCdfJetFinder.cxx:872
 AliCdfJetFinder.cxx:873
 AliCdfJetFinder.cxx:874
 AliCdfJetFinder.cxx:875
 AliCdfJetFinder.cxx:876
 AliCdfJetFinder.cxx:877
 AliCdfJetFinder.cxx:878
 AliCdfJetFinder.cxx:879
 AliCdfJetFinder.cxx:880
 AliCdfJetFinder.cxx:881
 AliCdfJetFinder.cxx:882
 AliCdfJetFinder.cxx:883
 AliCdfJetFinder.cxx:884
 AliCdfJetFinder.cxx:885
 AliCdfJetFinder.cxx:886
 AliCdfJetFinder.cxx:887
 AliCdfJetFinder.cxx:888
 AliCdfJetFinder.cxx:889
 AliCdfJetFinder.cxx:890
 AliCdfJetFinder.cxx:891
 AliCdfJetFinder.cxx:892
 AliCdfJetFinder.cxx:893
 AliCdfJetFinder.cxx:894
 AliCdfJetFinder.cxx:895
 AliCdfJetFinder.cxx:896
 AliCdfJetFinder.cxx:897
 AliCdfJetFinder.cxx:898
 AliCdfJetFinder.cxx:899
 AliCdfJetFinder.cxx:900
 AliCdfJetFinder.cxx:901
 AliCdfJetFinder.cxx:902
 AliCdfJetFinder.cxx:903
 AliCdfJetFinder.cxx:904
 AliCdfJetFinder.cxx:905
 AliCdfJetFinder.cxx:906
 AliCdfJetFinder.cxx:907