ROOT logo
//*-- Author : Jochen Markert 18.07.2007

#include "hparticletracksorter.h"
#include "hparticletool.h"
#include "hades.h"
#include "hruntimedb.h"
#include "hevent.h"
#include "hparticlecand.h"
#include "hemccluster.h"
#include "hparticledef.h"
#include "emcdef.h"
#include "htool.h"

#include <algorithm>
#include <iostream>
#include <iomanip>

using namespace std;
#define DEBUG 0

//_HADES_CLASS_DESCRIPTION
////////////////////////////////////////////////////////////////////////////////
//
//
// HParticleTrackSorter
//
// This Class loops over the HParticleCand category and
// analyze the objects. If the ECAL detector is used a loop over
// HEmcCluster is performed too to mark the clusters used by selected
// candidates. The Purpose of the procedure is to categorize
// the objects in the case of multiple usage of the single detector
// hits (RICH, InnerMDC, OuterMDC, META) as Double_t hit and to provide
// a hint which of the candidates is the best with respect to certain
// criteria (described by the bit flags below). If a HParticleCand
// object has been marked before by an other task with HParticleCand::kIsUsed it
// will be ignored by the sort process and stay unchanged. For the description
// of the differnt set/get functions of the HParticleCand objects see the
// documentation.
// To use HParticleTrackSorter without changing the flags of the HParticleCand objects
// (non persitent way inside a macro or reconstructor) one can call the
// backupFlags() function before and modification of the flags and
// restore the flags by calling restoreFlags() after the privat work has been
// done.
//
// For debugging purpose one cand set screen printouts to check the sorting
// (setDebug() + setPrintLevel(Int_t level)(level = 1-3)).
//
// The analysis of the HParticleCand objects is done in several steps:
//
// 1. Create the HParticleTrackSorter object and call the HParticleTrackSorter::Init()
//    To initialze the sorter object with category pointers etc ...
// 2. The internal structures are filled per event from the HParticleCand
//    category by calling the
//    Int_t HParticleTrackSorter::fill(Bool_t (*function)(HParticleCand* ))
//    function which expects a function pointer to a function taking an
//    HParticleCand pointer as argument and returning a Bool_t if the conditions on
//    the HParticleCand object are fullfilled (like for example in
//    Bool_t HParticleTrackSorter::selectHadrons(HParticleCand* pcand). For objects
//    failing the test the kIsRejected bit is set. The user can provide
//    his own selection functions to replace the build in selectLeptons(.....)
//    and selectHadrons(....) select functions.
//    The function pointer can be a pointer to a global function or member function
//    of an Object for example:
//
//    Bool_t select(HParticleCand* cand){             // global function
//             if ( put all selection criteria here ) return kTRUE;
//             else                                   return kFALSE;
//    } or
//
//    static Bool_t dummy::select(HParticleCand* cand){ // member function of object dummy
//                                                      // needs to be declared static !
//             if ( put all selection criteria here ) return kTRUE;
//             else                                   return kFALSE;
//    }
//    would be called in the code like
//
//    dummy d;
//    HParticleCand* p= new HParticleCand() // just to get an object
//    HParticleTrackSorter sorter;               // create the sorter object
//    // now we have 3 possibilities to call the
//    // test function
//    1. sorter.fill(select)        // global function
//    2. sorter.fill(dummy::select) // member function of object dummy (static call without object creation)
//    3. sorter.fill(d.select)      // member function of real object dummy
// 3. For all accepted objects the quality criteria and Double_t hit flags are
//    evaluated in within the selected sample and flagged to the HParticleCand objects.
//    The objects which fullfill the given minimum criteria the kIsAcceptedXXX bit is
//    flagged. For the RICH, InnerMDC, OuterMDC and META the kIsDoubleHitXXX bit are set
//    if the same detector hit index has been used more than once in the sample.
//    For the double hit rejection of Meta hits only the selected Meta hits
//    (ParticleCand::getSelectedMeta()) are taken into account. In the case of RPCClst+SHOWERHIT
//    The shower hit might be used by more than one candidate.
// 4. With the call of
//    Int_t HParticleTrackSorter::selectBest(Particle::ESwitch byQuality, Int_t byParticle)
//    the final selection on the objects is performed. The by the fill() procedure
//    selected objects are sorted by the selected quality criteria. Starting from the
//    best object one goes down the list and marks all objects with the kIsUsed and
//    kIsBestXXX bit if they do not reuse an already used detector hit index
//    (Double_t hit rejection).
//    For leptons the kIsLepton bit is flagged in addition. Individual detector hits
//    can be ignored in the selection procedure with the setIgnoreXXX(kTRUE) functions
//    (by default all hits are taken into account).
//    Particle::ESwitch byQuality can be used as quality criteria switch for
//
//          Particle::kIsBestHitRichSorter          (by number of pads)
//          Particle::kIsBestHitInnerMdcSorter      (by chi2, non fitted segments last)
//	    Particle::kIsBestHitOuterMdcSorter      (by chi2, non fitted segments last)
//          Particle::kIsBestHitMetaSorter          (by RK metamatch quality (or metaMatch quality if no RK , lower proirity))
//          Particle::kIsBestRKSorter               (by RK chi2, none fitted outer segments with lower priority)
//          Particle::kIsBestRKRKMETASorter         (by RK chi2 * RK META match quality)
//          Particle::kIsBestRKRKMETARadiusSorter   (by RK chi2 * RK META match radius)
//          Particle::kIsBestUserSorter             (by user provided sort function, no default function!)
//
//    user provide sort function:
//    static void   setUserSort(Bool_t    (*function)(candidateSort*, candidateSort*)
//    The user has to provide a function for sorting of the candidates:
//
//    Bool_t mysort(candidateSort* a, candidateSort* b) {
//      // this example does the same as kIsBestRKRKMETA
//      // sort in ascending order by RK_Chi2 * RK_META_match_Quality
//      // the function can use all data members of struct candidateSort defined in hparticletracksorter.h
//      return (a->RK_Chi2 * a->RK_META_match_Quality)  < (b->RK_Chi2 * b->RK_META_match_Quality);
//    }
//
//
// 5. The scheme allows for multiple selection following each other as
//    the already used objects are excluded from the next steps (for example first
//    select leptons and in the second iteration hadrons). The detector hit indicies
//    of the previous selection of kIsUsed marked objects are added to list of used
//    indicies by default not allowing to use a detector hit twice. To skip that
//    memory one can set setIgnorePreviousIndex(kTRUE).
// 6. After the analysis call HParticleTrackSorter::finalize(). After the selection
//    the user has to take into account only HParticleCand objects which  have
//    been flagged with kIsUsed. The objects slected by the Lepton selection will
//    be marked in addition with kIsLepton to identify them easily.
// 7. For lepton selections the internal selectLepton() function cand be manipulated.
//    With setRICHMatching(Particle::EMatch ,Float_t window) the matching between
//    inner MDC and RICH hit can be selected (default: HParticleTrackSorter::kUseRKRICHWindow,4 degree ).
// 8. For Better suppression of background the build in select functions can apply
//    a spacial match to the meta cell in y coordinate. setUseYMatching(kTRUE) (default)
//    makes use of a match in y in within 3 mm. ( see HParticleTool::isGoodMetaCell(HParticleCand* c,Float_t bound))
// 9. For Better suppression of background the build in select functions can require a valid beta measurement (beta>0)
//    setUseBeta(kTRUE) (default), for leptons in addition a low cut in beta can be set setBetaLetonCut(0.9) (default)
//10. For Better suppression of background the build in select functions can suppress fake measurement
//    setUseFakeRejection(kTRUE) (default)
//
//
//    A selection of unique lepton and hadron tracks could look like
//    the following example :
//
// void myMacro(Int_t nEvents,TString inputDir, TString inFile)
// {
//    Hades* myHades = new Hades;
//    gHades->setQuietMode(2);
//
//    HRootSource* source = new HRootSource();
//    source->setDirectory(((Text_t *)inputDir.Data()));
//    source->addFile((Text_t *)inFile.Data());
//
//    myHades->setDataSource(source);
//
//    if(!myHades->init()){
//	cout<<"Hades Init() failed!"<<endl;
//	exit(1);
//    }
//
//
//    //------------------------------------------------------------------------
//    // get category ointers and iterators
//    HCategory* pidCat  = (HCategory*)gHades->getCurrentEvent()->getCategory(catParticleCand);
//    HIterator* iter    = NULL;
//    if(pidCat)
//    {
//         iter = (HIterator *)pidCat->MakeIterator("native");
//    } else {
//        cout<<"Category Pointer is NULL!"<<endl;
//	exit(1);
//    }
//    //------------------------------------------------------------------------
//
//    --------------------------------------------------------------------------------------------
//    //--------------------------CONFIGURATION---------------------------------------------------
//    //At begin of the program (outside the event loop)
//    HParticleTrackSorter sorter;
//    //sorter.setDebug();                                            // for debug
//    //sorter.setPrintLevel(3);                                      // max prints
//    //sorter.setRICHMatching(Particle::kUseRKRICHWindowSorter,4.); // select matching RICH-MDC for selectLeptons() function
//    //sorter.setUseYMatching(kTRUE);                                // use meta cell build in select functions
//    //sorter.setMetaBoundary(4.0);                                   // match in 4 mm to meta cell
//    //sorter.setUseBeta(kTRUE);                                     // require beta >0 in build in select functions
//    //sorter.setTOFMAXCut( 60.);                                    // maximum tof [ns] allowd for match ( if beta is used)
//    //sorter.setBetaLeptonCut(0.9);                                 // lower beta cut for lepton select
//    //sorter.setUseFakeRejection(kTRUE);                            // reject fakes in  build in select functions (default kTRUE)
//    //sorter.setUseMETAQA(kTRUE);                                   // metaqa in  build in select functions (default kTRUE)
//    //sorter.setMETAQACut(3.0);                                     // metaqa cut for ahdrons+leptons
//    //sorter.setIgnoreInnerMDC();                                   // do not reject Double_t inner MDC hits
//    //sorter.setIgnoreOuterMDC();                                   // do not reject Double_t outer MDC hits
//    //sorter.setIgnoreMETA();                                       // do not reject Double_t META hits
//    //sorter.setIgnorePreviousIndex();                              // do not reject indices from previous selctions
//    sorter.init();                                                  // get catgegory pointers etc...
//    --------------------------------------------------------------------------------------------
//
//
//    //------------------------------------------------------------------------
//    // Event loop
//    for(Int_t i = 1;i < nEvents; i ++)
//    {
//        //----------break if last event is reached-------------
//        if(!gHades->eventLoop(1)) break;
//
//        //------------------------------------------------------------------------
//        // run your selection and set all flags
//
//        //------------------------------------------------------------------------
//        // clean vectors and index arrays
//        sorter.cleanUp();
//        //------------------------------------------------------------------------
//
//        // if you do not want to modify the flags
//        // for the following tasks etc (needs call of
//        // restoreFlags() at the end of the eventloop)
//        //sorter.backupFlags();
//
//        sorter.resetFlags(kTRUE,kTRUE,kTRUE,kTRUE);     // reset all flags for flags (0-28) ,reject,used,lepton
//        Int_t nCandLep     = sorter.fill(HParticleTrackSorter::selectLeptons);   // fill only good leptons
//        Int_t nCandLepBest = sorter.selectBest(Particle::kIsBestRKSorter,HParticleTrackSorter::kIsLeptonSorter);
//        Int_t nCandHad     = sorter.fill(HParticleTrackSorter::selectHadrons);   // fill only good hadrons (already marked good leptons will be skipped)
//        Int_t nCandHadBest = sorter.selectBest(Particle::kIsBestRKSorter,HParticleTrackSorter::kIsHadronSorter);
//        //------------------------------------------------------------------------
//
//        //------------------------------------------------------------------------
//        // analylize marked tracks (kIsUsed) ....
//    	  iter->Reset();
//        HParticleCand* pcand;
//        while((pcand = (HParticleCand*)iter->Next()) != NULL)
//	  {
//           // accept only objects marked as good
//	     if(!pcand->isFlagBit(HParticleCand::kIsUsed) ) continue;
//
//           // do other selection......
//        }
//        //------------------------------------------------------------------------
//
//        // if you do not want to modify the flags
//        // for the following tasks etc (needs call of
//        // backupFlags() at the beginning of the eventloop)
//        // sorter.restoreFlags();
//
//
//    }
//    //------------------------------------------------------------------------
//    //At the end of program
//    sorter.finalize();
// }
////////////////////////////////////////////////////////////////////////////////

ClassImp(sorter_setup)
ClassImp(HParticleTrackSorter)

Bool_t  HParticleTrackSorter::kDebug               = kFALSE;
Int_t   HParticleTrackSorter::printLevel           = 1;

Int_t   HParticleTrackSorter::kSwitchIndex         = 0;
Int_t   HParticleTrackSorter::kSwitchQuality       = 0;
Int_t   HParticleTrackSorter::kSwitchParticle      = 0;
Int_t   HParticleTrackSorter::kSwitchRICHMatching  = Particle::kUseRKRICHWindowSorter;
Float_t HParticleTrackSorter::fRICHMDCWindow       = 4.;

Bool_t  HParticleTrackSorter::kIsFieldOn           = kTRUE;
Bool_t  HParticleTrackSorter::kIgnoreRICH          = kFALSE;
Bool_t  HParticleTrackSorter::kIgnoreInnerMDC      = kFALSE;
Bool_t  HParticleTrackSorter::kIgnoreOuterMDC      = kFALSE;
Bool_t  HParticleTrackSorter::kIgnoreMETA          = kFALSE;
Bool_t  HParticleTrackSorter::kIgnorePreviousIndex = kFALSE;
Bool_t  HParticleTrackSorter::kUseYMatching        = kTRUE;
Bool_t  HParticleTrackSorter::kUseYMatchingScaling = kTRUE;
Float_t HParticleTrackSorter::fMetaBoundary        = 4.0;
Bool_t  HParticleTrackSorter::kUseBeta             = kTRUE;
Float_t HParticleTrackSorter::fBetaLepCut          = 0.9;
Bool_t  HParticleTrackSorter::kUseFakeRejection    = kTRUE;
Bool_t  HParticleTrackSorter::kUseMETAQA           = kTRUE;
Float_t HParticleTrackSorter::fMETAQACut           = 3.0;
Float_t HParticleTrackSorter::fTOFMAXCut           = 60.0;
Float_t HParticleTrackSorter::fDeltaSegCut         = 2.0;

Bool_t (*HParticleTrackSorter::pUserSort)(candidateSort*, candidateSort*) = NULL;

HParticleTrackSorter::HParticleTrackSorter(void)
    : TNamed("ParticleSorter","PID track candidate sorter")
{
    clear();
    setup.setDefault();
}

HParticleTrackSorter::HParticleTrackSorter(TString name, TString title)
: TNamed(name, title)
{
    clear();
    setup.setDefault();
}

HParticleTrackSorter::~HParticleTrackSorter(void)
{
    if(nameIndex)           delete [] nameIndex;
    if(nameQuality)         delete [] nameQuality;
    if(iterParticleCandCat) delete iterParticleCandCat;
}

void HParticleTrackSorter::clear()
{
    // initialize the member variables
    nt                   = NULL;
    isSimulation         = kFALSE;
    nameIndex            = 0;
    nameQuality          = 0;
    pParticleCandCat     =  NULL;
    pEmcClusterCat       =  NULL;
    iterParticleCandCat  =  NULL;
    fill_Iteration       =  0;
    selectBest_Iteration =  0;
    currentEvent         = 100000000;
    fEvent               = NULL;
}

Bool_t HParticleTrackSorter::init(HRecEvent* evt)
{
    // get the pointer to the HParticleCand categegory and create the
    // corresponding iterator. HRecEvent* evt == 0 (default) the current event
    // will be used otherise the event provided.

    if(evt) { fEvent = evt; }
    else    { fEvent = (HRecEvent*)gHades->getCurrentEvent(); }

    pParticleCandCat = fEvent->getCategory(catParticleCand);
    if (!pParticleCandCat) {
	Error("HParticleTrackSorter::init()"," No HParticleCand Input -> Switching HParticleTrackSorter OFF");
    }
    else {
        if(iterParticleCandCat) delete iterParticleCandCat;
	iterParticleCandCat = (HIterator *) pParticleCandCat->MakeIterator();
    }
    pEmcClusterCat   = fEvent->getCategory(catEmcCluster);

    HCategory* catKine = fEvent->getCategory(catGeantKine);
    if(catKine) isSimulation = kTRUE;

    if(kDebug) {

	Int_t size;
	nameIndex   = HTool::parseString("kIsIndexRICH,kIsIndexInnerMDC,kIsIndexOuterMDC,"
					 "kIsIndexTOF,kIsIndexSHOWER,kIsIndexRPC,kIsIndexMETA",size,",",kFALSE);
        cout<<"-------------------------------------------"<<endl;
	cout<<"HParticleTrackSorter::init() :"<<endl;
	cout<<"index :: " <<endl;
	for(Int_t i = 0; i < size; i ++){
	    cout<<i <<" "<<nameIndex[i].Data()<<endl;
	}
	nameQuality = HTool::parseString("kIsBestHitRICH,kIsBestHitInnerMDC,kIsBestHitOuterMDC,kIsBestHitMETA,"
					 "kIsBestRK,kIsBestRKRKMETA,kIsBestRKRKMETARadius",size,",",kFALSE);
        for(Int_t i = 0; i < size; i ++){
	    cout<<i + SORTER_NUMBER_OF_INDEX<<" "<<nameQuality[i].Data()<<endl;
	}
	cout<<"-------------------------------------------"<<endl;
    }

     return kTRUE;
}

Bool_t HParticleTrackSorter::finalize(void)
{
    // clean up the memory and
    // write ntuple to file (if it exist)
    cleanUp(kTRUE);
    return kTRUE ;
}

Bool_t  HParticleTrackSorter::isDeltaSegOk(HParticleCand* cand)
{
    // returns kTRUE if inner and outer segment
    // agree in theta, phi inside fDeltaSegCut

    if(cand->getTheta2() ==-1) return kFALSE;

    if((fabs(cand->getPhi()   - cand->getPhi2())   < fDeltaSegCut) &&
       (fabs(cand->getTheta() - cand->getTheta2()) < fDeltaSegCut)
      ) return kTRUE;

    return kFALSE;

}


Int_t HParticleTrackSorter::clearVector( vector<candidateSort*>& all_candidateSorts)
{
    // delete all candidate structs in vector
    // from heap and call the clear function of
    // vector afterwards.
    for(UInt_t i = 0; i < all_candidates.size() ; ++ i)
    {
	candidateSort* cand = all_candidates[i];
        delete cand;
    }
    all_candidates.clear();
    return 0;
}


Bool_t HParticleTrackSorter::cmpIndex(candidateSort* a, candidateSort* b)
{
    // sorts the candidate list ascending by the index
    // of the detector hit of interest.

    switch (kSwitchIndex) {
    case Particle::kIsIndexRICHSorter:
	// sort in ascending order by Rich index
	return a->ind_RICH     < b->ind_RICH;
    case Particle::kIsIndexInnerMDCSorter:
	// sort in ascending order by inner MDC segment index
	return a->ind_innerMDC < b->ind_innerMDC;
    case Particle::kIsIndexOuterMDCSorter:
	// sort in ascending order by outer MDC segment index
	return a->ind_outerMDC < b->ind_outerMDC;
    case Particle::kIsIndexTOFSorter:
	// sort in ascending order by TOF index
	return a->ind_TOF      < b->ind_TOF;
    case Particle::kIsIndexSHOWERSorter:
	// sort in ascending order by SHOWER index
	return a->ind_SHOWER   < b->ind_SHOWER;
    case Particle::kIsIndexRPCSorter:
	// sort in ascending order by RPC index
	return a->ind_RPC      < b->ind_RPC;

    default:
	::Error("HParticleTrackSorter::CmpIndex()","Unknown switch statement %i!",kSwitchIndex);
	return kFALSE;

    }
}

Bool_t HParticleTrackSorter::cmpQuality(candidateSort* a, candidateSort* b)
{
    // sorts the candidate list by the quality critera
    // of interest (ascending or decending depending on
    // the quantity). The best object is always the first.

    switch (kSwitchQuality) {
    case Particle::kIsBestHitRICHSorter:
	// sort in ascending order by inner MDC chi2
	return a->nPadsRich             < b->nPadsRich;
    case Particle::kIsBestHitInnerMDCSorter:
	// sort in ascending order by inner MDC chi2
	return a->innerMDCChi2          < b->innerMDCChi2;
    case Particle::kIsBestHitOuterMDCSorter:
	// sort in ascending order by outer MDC chi2
	return a->outerMDCChi2          < b->outerMDCChi2;
    case Particle::kIsBestHitMETASorter:
	// sort in ascending order by RK METAMATCH QUALITY
	return a->RK_META_match_Quality < b->RK_META_match_Quality;
    case Particle::kIsBestRKSorter:
        // sort in ascending order by RK chi2
	return a->RK_Chi2               < b->RK_Chi2;
    case Particle::kIsBestRKRKMETASorter:
	// sort in ascending order by RK_Chi2 * RK_META_match_Quality
        return (a->RK_Chi2 * a->RK_META_match_Quality)  < (b->RK_Chi2 * b->RK_META_match_Quality);
    case Particle::kIsBestRKRKMETARadiusSorter:
	// sort in ascending order by RK_Chi2 * RK_META_match_Quality
        return (a->RK_Chi2 * a->RK_META_match_Radius)  < (b->RK_Chi2 * b->RK_META_match_Radius);
    case Particle::kIsBestUserSorter:
	// sort in ascending order by user function
	if(pUserSort == NULL) {
            ::Error("HParticleTrackSorter::CmpQuality()","kIsBestUser is set but user function pointer i NULL!");
	    return kTRUE;
	}
        else return (*pUserSort)(a,b);
    default:
	::Error("HParticleTrackSorter::CmpQuality()","Unknown switch statement %i!",kSwitchQuality);
        return kFALSE;

    }
}


Bool_t HParticleTrackSorter::rejectIndex(candidateSort* a, Particle::ESwitch byIndex, Int_t& index)
{
    // if function return kTRUE the Object will be rejected
    // rejects candidates from the actual list, where the
    // detector hit of interest is not defined.

    switch (byIndex) {
    case Particle::kIsIndexRICHSorter:
	// NO RICH
	return ( index = a->ind_RICH )     < 0;
    case Particle::kIsIndexInnerMDCSorter:
	// NO inner MDC segment
	return ( index = a->ind_innerMDC ) < 0;
    case Particle::kIsIndexOuterMDCSorter:
	// NO outer MDC segment
	return ( index = a->ind_outerMDC ) < 0;
    case Particle::kIsIndexTOFSorter:
	// NO TOF
	return ( index = a->ind_TOF )      < 0;
    case Particle::kIsIndexSHOWERSorter:
	// NO SHOWER
	return ( index = a->ind_SHOWER )   < 0;
    case Particle::kIsIndexRPCSorter:
	// NO RPC
	return ( index = a->ind_RPC )      < 0;
    default:
	::Error("HParticleTrackSorter::rejectIndex()","Unknown switch statement %i!",byIndex);
        return kFALSE;
    }
}
Bool_t HParticleTrackSorter::rejectQuality(candidateSort* a, Particle::ESwitch byQuality)
{
    // if function return kTRUE the Object will be rejected
    // rejects the candidates from the actual list if the quality
    // of interest does not match the mimnum requirement.

    switch (byQuality) {
    case Particle::kIsBestHitRICHSorter:
	// RICH ring
	return a->ind_RICH < 0;
    case Particle::kIsBestHitInnerMDCSorter:
	// NO fitted inner MDC segments
	return a->innerMDCChi2 < 0;
    case Particle::kIsBestHitOuterMDCSorter:
	// NO fitted outer MDC segments
	return (a->outerMDCChi2 < 0 || a->ind_outerMDC < 0);
    case Particle::kIsBestHitMETASorter:
	// NO META
	return (a->ind_TOF < 0 && a->ind_SHOWER < 0 && a->ind_RPC < 0);
    case Particle::kIsBestRKSorter:
	// NO outer MDC segment or NO RK
	return (a->RK_Chi2 < 0);
    case Particle::kIsBestRKRKMETASorter:
	// NO META or NO RK
        return (a->ind_TOF < 0 && a->ind_SHOWER < 0 && a->ind_RPC < 0) || (a->RK_Chi2 < 0);
    case Particle::kIsBestRKRKMETARadiusSorter:
	// NO META or NO RK
        return (a->ind_TOF < 0 && a->ind_SHOWER < 0 && a->ind_RPC < 0) || (a->RK_Chi2 < 0);
    case Particle::kIsBestUserSorter:
	// NO META or NO RK
        return (a->ind_TOF < 0 && a->ind_SHOWER < 0 && a->ind_RPC < 0) || (a->RK_Chi2 < 0);

    default:
	::Error("HParticleTrackSorter::rejectQuality()","Unknown switch statement %i!",byQuality);
        return kFALSE;

    }
}

Int_t HParticleTrackSorter::flagAccepted(vector <candidateSort*>& all_candidates, Particle::ESwitch byQuality)
{
    // flag all candidates which does satisfy the quality
    // criteria from input with kIsAcceptedXXX bit. If kDebug
    // is set and printLevel >=2 the statistics for the accepted
    // objects is printed.


    Int_t ctAll           = 0;
    Int_t ctAcceptQuality = 0;
    HParticleCand* pcand;
    for(UInt_t i = 0; i < all_candidates.size(); ++ i)
    {
        ++ ctAll;
	candidateSort* cand = all_candidates[i];
	if(rejectQuality(cand, byQuality))      continue;
	++ ctAcceptQuality;

	pcand = dynamic_cast<HParticleCand*>(pParticleCandCat->getObject(cand->ind_Cand));
	if(!pcand) {
	    Error("HParticleTrackSorter::flagAccepted()","NULL pointer retrieved for HParticleCand at index %i",cand->ind_Cand);
	    continue;
	}
	switch (byQuality) {
	case Particle::kIsBestHitRICHSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedHitRICH);
            //------------------------------------------------
            // check the ring correlation
	    if(pcand->getChi2() >= 0){
		if( pcand->getRingCorr() == Particle::kIsRICHRK ||
		    pcand->getRingCorr() == (Particle::kIsRICHRK|Particle::kIsRICHMDC)
		  ) pcand->setFlagBit(Particle::kIsAcceptedRKRICH);
	    }
	    if( pcand->getRingCorr() == Particle::kIsRICHMDC ||
	        pcand->getRingCorr() == (Particle::kIsRICHRK|Particle::kIsRICHMDC)
	      ) pcand->setFlagBit(Particle::kIsAcceptedHitRICHMDC);
            //------------------------------------------------
	    break;
	case Particle::kIsBestHitInnerMDCSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedHitInnerMDC);
	    break;
	case Particle::kIsBestHitOuterMDCSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedHitOuterMDC);
	    break;
	case Particle::kIsBestHitMETASorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedHitMETA);
            //------------------------------------------------
            // check if it was accepted for RK
	    if(pcand->getChi2() >= 0) {
		pcand-> setFlagBit(Particle::kIsAcceptedRKMETA);
	    }
            //------------------------------------------------
	    break;
	case Particle::kIsBestRKSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedRK);
	    break;
	case Particle::kIsBestRKRKMETASorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedRKRKMETA);
	    break;
	case Particle::kIsBestRKRKMETARadiusSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedRKRKMETARadius);
	    break;
	case Particle::kIsBestUserSorter:
	    pcand -> setFlagBit(Particle::kIsAcceptedUser);
	    break;
	default :

	    Error("HParticleTrackSorter::flagAccepted()","Unknown switch statement : %i!",byQuality);
	    exit(1);
	    break;
	}

    }
    if(kDebug && printLevel >= 2){
	cout<<"flagAccepted(): "
	    <<"s: "                  <<selectBest_Iteration
	    <<", f: "                <<fill_Iteration
	    <<", kSwitchParticle: "  <<kSwitchParticle
	    <<", quality "           <<nameQuality[byQuality-SORTER_NUMBER_OF_INDEX].Data()
	    <<", candidate= "        <<ctAll
	    <<", acc. quality = "    <<ctAcceptQuality
	    <<endl;

    }
    return 0;
}

Int_t HParticleTrackSorter::flagDouble(vector <candidateSort*>& all_candidates, Particle::ESwitch byIndex)
{
    // sort vector in ascending order by the particular index byindex. Loops
    // over the vector and extract objects with same detetector index and marks them
    // as Double_t hits in HParticleCand (kIsDoubleHitXXX bit). Non valid detecor
    // hits (index -1) are not taken into account. If kDebug
    // is set and printLevel >=2 the statistics for the accepted
    // objects is printed.

    kSwitchIndex   = byIndex;
    sort(all_candidates.begin(), all_candidates.end(), this->cmpIndex);


    Int_t index;
    Int_t old_index =-99;
    vector <candidateSort*> daughter;

    Int_t ctAll          = 0;
    Int_t ctAcceptIndex  = 0;
    Int_t ctAcceptDoubles= 0;
    Int_t ctCheck        = 0;

    for(UInt_t i = 0; i < all_candidates.size(); ++ i)
    {
        ++ ctAll;
	candidateSort* cand = all_candidates[i];
	if(rejectIndex  (cand, byIndex, index)) continue;
        ++ ctAcceptIndex;
	if(index == old_index || old_index == -99)
	{ // first time  or  still same index
	    daughter.push_back(cand);
	} else {
            // finalize list
	    if(!daughter.empty())
	    {
                setFlagsDouble(daughter,byIndex);
		if(daughter.size() >= 2){
		    ctAcceptDoubles += daughter.size();
		}
		ctCheck += daughter.size();
		daughter.clear();
	    }
	    daughter.push_back(cand);
	}
	old_index = index;
    }
    // finalize list last group
    if(!daughter.empty())
    {
	setFlagsDouble(daughter, byIndex);
	if(daughter.size() >= 2){
	    ctAcceptDoubles += daughter.size();
	}
	ctCheck += daughter.size();
    }
    if(kDebug && printLevel >= 2){
	cout<<"flagDouble(): "
	    <<"s: "                <<selectBest_Iteration
	    <<", f: "              <<fill_Iteration
	    <<", kSwitchParticle: "<<kSwitchParticle
	    <<", index "           <<nameIndex  [byIndex].Data()
	    <<", candidate= "      <<ctAll
	    <<", acc. index = "    <<ctAcceptIndex
            <<", acc. doubles = "  <<ctAcceptDoubles
            <<", check = "         <<ctCheck
	    <<endl;

    }
    return 0;

}

Int_t HParticleTrackSorter::setFlagsDouble(vector<candidateSort*>& daughter, Particle::ESwitch byIndex)
{
    // flag the HParticleCand objects with the rank of the Double_t hit criteria as output
    // form the sorting algorithms. Expects as input a vector of candidates with the
    // same hit index.

    if(daughter.empty()   ) return 0;    // nothing to do
    if(daughter.size() < 2) return 0;    // no Double_t hit


    HParticleCand* cand;
    for(UInt_t i = 0; i < daughter.size(); ++ i)
    {  // loop over list and set flags
	cand = dynamic_cast<HParticleCand*>(pParticleCandCat->getObject(daughter[i]->ind_Cand));
	if(!cand) {
	    Error("HParticleTrackSorter::setFlagsDouble()","NULL pointer retrieved for HParticleCand at index %i",daughter[i]->ind_Cand);
	    continue;
	}
	switch (byIndex)
	{
	case Particle::kIsIndexRICHSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitRICH);
	    break;
	case Particle::kIsIndexInnerMDCSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitInnerMDC);
	    break;
	case Particle::kIsIndexOuterMDCSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitOuterMDC);
	    break;
	case Particle::kIsIndexTOFSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitMETA);
	    break;
	case Particle::kIsIndexSHOWERSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitMETA);
	    break;
	case Particle::kIsIndexRPCSorter:
	    cand -> setFlagBit(Particle::kIsDoubleHitMETA);
	    break;
	default :

	    Error("HParticleTrackSorter::setFlagsDouble()","Unknown switch statement : %i!",byIndex);
	    exit(1);
	    break;

	}
    }

    return 0;
}

Bool_t HParticleTrackSorter::flagEmcClusters()
{
    // loop over EmcCluster cat and reset IsUsedInParticleCand.
    // loop over all candidates and flag emc clusters used in
    // selected candidates only
    if(pEmcClusterCat&&pParticleCandCat)
    {
        HEmcCluster* emc = 0;
	Int_t n = pEmcClusterCat->getEntries();
	for(Int_t i = 0; i < n; i ++){
            emc = (HEmcCluster*) pEmcClusterCat->getObject(i);
            if(emc) emc -> unsetIsUsedInParticleCand();
	}

        HParticleCand* cand = 0;
	n = pParticleCandCat->getEntries();
	for(Int_t i = 0; i < n; i ++){
	    cand = dynamic_cast<HParticleCand*>(pParticleCandCat->getObject(i));
	    if(cand->isFlagBit(Particle::kIsUsed) ){
		Int_t ind = cand->getEmcInd();
		if(ind>=0){
		    emc = (HEmcCluster*) pEmcClusterCat->getObject(ind);
                    if(emc)  emc->setIsUsedInParticleCand();
		}
	    }
	}
    }
     return kTRUE;
}


Int_t HParticleTrackSorter::fillInput(vector < candidateSort* >& all_candidates)
{
    // loop over HParticleCand category and fill a vector of temporary objects
    // which will be used for sorting and flagging of the objects in the end.
    // The candidate objects are initialzied with values of interest in a well
    // defined way, which makes the sorting and rejecting of objects later on easier.

    HParticleCand* pCand ;

    candidateSort* cand;
    iterParticleCandCat->Reset();
    while ((pCand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0 )
    {

	if(pCand->isFlagBit(Particle::kIsUsed))
	{
	    if(!kIgnorePreviousIndex && fill_Iteration == 2)
	    {
		//-------------------------------------------------------------------
		// add detector indices of already used objects to
		// the current index list
		if(pCand->getRichInd()     >= 0 &&
		   (find(index_RICH.begin(),index_RICH.end()        ,pCand->getRichInd())     == index_RICH.end()) ){
		    index_RICH.push_back(pCand->getRichInd());
		}
		if(pCand->getInnerSegInd() >= 0 &&
		   (find(index_InnerMDC.begin(),index_InnerMDC.end(),pCand->getInnerSegInd()) == index_InnerMDC.end()) ){
		    index_InnerMDC.push_back(pCand->getInnerSegInd());
		}
		if(pCand->getOuterSegInd() >= 0 &&
		   (find(index_OuterMDC.begin(),index_OuterMDC.end(),pCand->getOuterSegInd()) == index_OuterMDC.end()) ){
		    index_OuterMDC.push_back(pCand->getOuterSegInd());
		}

		//-------------------------------------------------------------------
                // TOF Clst + TOF hit special case
		if(pCand->getSelectedMeta() == Particle::kTofClst && pCand->getTofClstInd()      >= 0 &&
		   (find(index_TOFClst.begin(),index_TOFClst.end()  ,pCand->getTofClstInd())  == index_TOFClst.end()) ){
		    index_TOFClst.push_back(pCand->getTofClstInd());
		}
		if((pCand->getSelectedMeta() == Particle::kTofHit1 || pCand->getSelectedMeta() == Particle::kTofHit2) &&
		   pCand->getTofHitInd()       >= 0 &&
		   (find(index_TOFHit.begin(),index_TOFHit.end()    ,pCand->getTofHitInd()+1000)  == index_TOFHit.end()) ){
		    index_TOFHit.push_back(pCand->getTofHitInd()+1000);
		}
		//-------------------------------------------------------------------
		if(pCand->getSelectedMeta() == Particle::kShowerHit && pCand->getShowerInd()   >= 0 &&
		   (find(index_SHOWER.begin(),index_SHOWER.end()    ,pCand->getShowerInd())   == index_SHOWER.end()) ){
		    index_SHOWER.push_back(pCand->getShowerInd());
		}

		if(pCand->getSelectedMeta() == Particle::kRpcClst && pCand->getRpcInd()      >= 0 &&
		   (find(index_RPC.begin(),index_RPC.end()          ,pCand->getRpcInd())      == index_RPC.end()) ){
		    index_RPC.push_back(pCand->getRpcInd());
		}

	    }
	    //-------------------------------------------------------------------

	    continue; // skip objects which are marked
	}
	if(pCand->isFlagBit(Particle::kIsRejected))  continue; // skip objects which are marked

	for(Int_t i = 0; i < 29; i ++) { pCand->unsetFlagBit(i);  }  // clean condition, only reject/used/lepton bit stays!

        //-------------------------------------------------------------------
	// fill the structur objects with needed
        // information

	cand = new candidateSort;
        cand->cand                  = pCand;
	cand->ind_Cand              = pCand->getIndex();

	cand->ind_RICH              = pCand->getRichInd();
	cand->ind_innerMDC          = pCand->getInnerSegInd();
	cand->ind_outerMDC          = pCand->getOuterSegInd();
	cand->ind_TOF               = (pCand->getTofClstInd() >= 0)? pCand->getTofClstInd() : (pCand->getTofHitInd() >= 0)? pCand->getTofHitInd() + 1000 : -1 ; // special case tof
	cand->ind_SHOWER            = pCand->getShowerInd();
	cand->ind_RPC               = pCand->getRpcInd();
	cand->selectedMeta          = pCand->getSelectedMeta();
        cand->nPadsRich             = pCand->getRingNumPads();
	cand->innerMDCChi2          = pCand->getInnerSegmentChi2();
	cand->outerMDCChi2          = pCand->getOuterSegmentChi2();
                                      
	cand->sec                   =  pCand->getSector();
	cand->beta                  =  pCand->getBeta();
	cand->mom                   =  pCand->getMomentum();
	cand->MDC_dEdx              =  pCand->getMdcdEdx();
	cand->nMdcLayerInnerSeg     =  pCand->getNLayer(0);
        cand->nMdcLayerOuterSeg     =  pCand->getNLayer(1);
	cand->RK_Chi2               = (pCand->getOuterSegmentChi2() < 0 && pCand->getChi2() >= 0 )? pCand->getChi2() + 100 : pCand->getChi2(); // lower priority sorting
        cand->RICH_RK_Corr          =  pCand->getRingCorr();
	cand->RK_RICH_match_Quality = (pCand->getRichInd() >= 0)?  pCand->getRichMatchingQuality(): -1;
	cand->RK_META_match_Quality = (finite(pCand->getMetaMatchQuality()))? pCand->getMetaMatchQuality() : -1;
	cand->RK_META_match_Radius  = (finite(pCand->getMetaMatchRadius()))? pCand->getMetaMatchRadius() : -1;
	cand->RK_META_dx            = (finite(pCand->getRkMetaDx()))? pCand->getRkMetaDx() : -1;
	cand->RK_META_dy            = (finite(pCand->getRkMetaDy()))? pCand->getRkMetaDy() : -1;
	cand->RK_META_match_Shr_Quality = pCand->getMetaMatchQualityShower();
        cand->innerMDCTheta = pCand->getTheta();
        cand->innerMDCPhi   = pCand->getPhi();
        cand->outerMDCTheta = pCand->getTheta2();
        cand->outerMDCPhi   = pCand->getPhi2();

	cand->userSort1=0;
	cand->userSort2=0;
	cand->userSort3=0;


	if(pCand->getSelectedMeta() == kTofClst ||
	   pCand->getSelectedMeta() == kTofHit1 ||
	   pCand->getSelectedMeta() == kRpcClst
	  ) {
	    cand->MetaModule = pCand->getMetaModule(0);
	    cand->MetaCell   = pCand->getMetaCell(0);
	} else if (pCand->getSelectedMeta() == kTofHit2){
	    cand->MetaModule = pCand->getMetaModule(1);
	    cand->MetaCell   = pCand->getMetaCell(1);
	} else {
	    cand->MetaModule = -1;
            cand->MetaCell   = -1;
	}

	if(cand->RK_Chi2 < 0){
	    cand->RK_META_match_Quality += 100; // lower priority in sorting
	    cand->RK_META_match_Radius  += 1000; // lower priority in sorting
	}

	//-------------------------------------------------------------------


	all_candidates.push_back(cand);
    }
    return 0;
}

void HParticleTrackSorter::printCand(candidateSort* cand, Int_t i,TString spacer)
{
    // print one candidate object and the flags which have been
    // already set to HParticleCand
    cout<<spacer.Data()<<i<<" --------------------------------------------"<<endl;
    cout<<spacer.Data()<<"cand " <<cand->ind_Cand
	<<", R: "  <<cand->ind_RICH
	<<", iM: " <<cand->ind_innerMDC
	<<", oM: " <<cand->ind_outerMDC
	<<", T: "  <<cand->ind_TOF
	<<", S: "  <<cand->ind_SHOWER
	<<", RPC: "<<cand->ind_RPC
	<<", selMeta: "<<cand->selectedMeta
	<<endl;
    cout<<spacer.Data()
	<<"Npads:"    <<cand->nPadsRich
	<<", ichi2: " <<cand->innerMDCChi2
	<<", ochi2: " <<cand->outerMDCChi2
	<<", RKM: "   <<cand->RK_META_match_Quality
	<<", RKMSHR: "<<cand->RK_META_match_Shr_Quality
	<<", RK: "    <<cand->RK_Chi2
	<<", RKRKMETA: "  <<(cand->RK_Chi2 != -1 && cand->RK_META_match_Quality != -1 ? cand->RK_Chi2*cand->RK_META_match_Quality : -1)
	<<", RKRKMETARadius: "  <<(cand->RK_Chi2 != -1 && cand->RK_META_match_Radius != -1 ? cand->RK_Chi2*cand->RK_META_match_Radius : -1)
	<<", isRingRK: "<<(Int_t) cand->RICH_RK_Corr
	<<endl;
    HParticleCand* pcand = dynamic_cast<HParticleCand*>(pParticleCandCat->getObject(cand->ind_Cand));
    if(cand)
    {
	cout<<spacer.Data()<<flush;
	pcand->printFlags();
    } else {
	::Error("HParticleTrackSorter::printCand()","NULL pointer retrieved for HParticleCand Index %i",cand->ind_Cand);

    }

}
void HParticleTrackSorter::printEvent(TString comment)
{
    Int_t ctcand = 0;

    cout<<"----------------------------------------"<<endl;
    cout<<"Event kSwitchParticle: "<<kSwitchParticle<<" "
 	<<setw(6)                  <<fEvent->getHeader()->getEventSeqNumber()
	<<"s: "                    <<selectBest_Iteration
	<<", f: "                  <<fill_Iteration
	<<" "                      <<comment.Data()
	<<endl;
    HParticleCand* pcand;
    iterParticleCandCat->Reset();
    while ((pcand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0 )
    {
	cout <<setw(6)<<ctcand<<" "<<flush;
	pcand->print();
	++ ctcand;
    }
}
void HParticleTrackSorter::printSetup()
{
    // print seeting of this sorter:
    // column 1 : backup
    // column 2 : setup
    // column 3 : actual

    cout<<"--------------------------------------------------------------------"<<endl;
    cout<<"HParticleTrackSorter: evt "<<gHades->getEventCounter()<<endl;
    cout<<"                        backup  setup  actual"<<endl;
    cout<<"kSwitchRICHMatching  "<<setup_backup.kSwitchRICHMatching <<" : "<<setup.kSwitchRICHMatching <<" : "<<kSwitchRICHMatching <<endl;
    cout<<"fRICHMDCWindow       "<<setup_backup.fRICHMDCWindow      <<" : "<<setup.fRICHMDCWindow      <<" : "<<fRICHMDCWindow      <<endl;
    cout<<"kIsFieldOn           "<<setup_backup.kIsFieldOn          <<" : "<<setup.kIsFieldOn          <<" : "<<kIsFieldOn          <<endl;
    cout<<"kIgnoreRICH          "<<setup_backup.kIgnoreRICH         <<" : "<<setup.kIgnoreRICH         <<" : "<<kIgnoreRICH         <<endl;
    cout<<"kIgnoreInnerMDC      "<<setup_backup.kIgnoreInnerMDC     <<" : "<<setup.kIgnoreInnerMDC     <<" : "<<kIgnoreInnerMDC     <<endl;
    cout<<"kIgnoreOuterMDCB     "<<setup_backup.kIgnoreOuterMDC     <<" : "<<setup.kIgnoreOuterMDC     <<" : "<<kIgnoreOuterMDC     <<endl;
    cout<<"kIgnoreMETA          "<<setup_backup.kIgnoreMETA         <<" : "<<setup.kIgnoreMETA         <<" : "<<kIgnoreMETA         <<endl;
    cout<<"kIgnorePreviousIndex "<<setup_backup.kIgnorePreviousIndex<<" : "<<setup.kIgnorePreviousIndex<<" : "<<kIgnorePreviousIndex<<endl;
    cout<<"pUserSort            "<<setup_backup.pUserSort           <<" : "<<setup.pUserSort           <<" : "<<pUserSort           <<endl;
    cout<<"kUseYMatching        "<<setup_backup.kUseYMatching       <<" : "<<setup.kUseYMatching       <<" : "<<kUseYMatching       <<endl;
    cout<<"kUseYMatchingScaling "<<setup_backup.kUseYMatchingScaling<<" : "<<setup.kUseYMatchingScaling<<" : "<<kUseYMatchingScaling<<endl;
    cout<<"fMetaBoundary        "<<setup_backup.fMetaBoundary       <<" : "<<setup.fMetaBoundary       <<" : "<<fMetaBoundary       <<endl;
    cout<<"kUseBeta             "<<setup_backup.kUseBeta            <<" : "<<setup.kUseBeta            <<" : "<<kUseBeta            <<endl;
    cout<<"fBetaLepCut          "<<setup_backup.fBetaLepCut         <<" : "<<setup.fBetaLepCut         <<" : "<<fBetaLepCut         <<endl;
    cout<<"kUseFakeRejection    "<<setup_backup.kUseFakeRejection   <<" : "<<setup.kUseFakeRejection   <<" : "<<kUseFakeRejection   <<endl;
    cout<<"kUseMETAQA           "<<setup_backup.kUseMETAQA          <<" : "<<setup.kUseMETAQA          <<" : "<<kUseMETAQA          <<endl;
    cout<<"fMetaBoundary        "<<setup_backup.fMetaBoundary       <<" : "<<setup.fMetaBoundary       <<" : "<<fMetaBoundary       <<endl;
    cout<<"fMETAQACut           "<<setup_backup.fMETAQACut          <<" : "<<setup.fMETAQACut          <<" : "<<fMETAQACut          <<endl;
    cout<<"fTOFMAXCut           "<<setup_backup.fTOFMAXCut          <<" : "<<setup.fTOFMAXCut          <<" : "<<fTOFMAXCut          <<endl;
    cout<<"fDeltaSegCut         "<<setup_backup.fDeltaSegCut        <<" : "<<setup.fDeltaSegCut        <<" : "<<fDeltaSegCut        <<endl;
    cout<<"--------------------------------------------------------------------"<<endl;

}

void HParticleTrackSorter::resetFlags(Bool_t flag, Bool_t reject, Bool_t used, Bool_t lepton)
{
    // reset all flags of HParticleCand by looping over the
    // HParticleCand category. The flags kIsRejected, kIsUsed
    // and kIsLepton are special for selection and handles independend
    // from the other flags. The lower flags of objects marked
    // with kIsUsed cand be only reste if used == kRTUE.

    if(!iterParticleCandCat) return;
    HParticleCand* pCand;
    iterParticleCandCat->Reset();
    while ((pCand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0) // begin of ParticleCand iterator
    {
	if(flag){
	    if((pCand->isFlagBit(Particle::kIsUsed) && used)  ||
	       !pCand->isFlagBit(Particle::kIsUsed) ){
		// we want to keep the flags from previous selection
		// clean condition, only reject,used and lepton bit stays!
		for(Int_t i=0;i<Particle::kIsLepton;i++) { pCand->unsetFlagBit(i);  }  
	    }
	}
	if(reject) pCand->unsetFlagBit(Particle::kIsRejected);
	if(used)   pCand->unsetFlagBit(Particle::kIsUsed);
	if(lepton) pCand->unsetFlagBit(Particle::kIsLepton);
    }
}
void HParticleTrackSorter::selection(Bool_t (*function)(HParticleCand* ))
{
    // loops over the HParticleCand category and marks all objects
    // which does not fullfill the selection criteria with the reject flag.
    // Objects marked as kIsUsed are skipped.
    if(!iterParticleCandCat) return;
    HParticleCand* pCand;
    iterParticleCandCat->Reset();
    while ((pCand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0) // begin of ParticleCand iterator
    {
	if( pCand->isFlagBit(Particle::kIsUsed)) { continue;} // do not touch already used objects!
	if(!pCand->select(function)) {pCand->setFlagBit(Particle::kIsRejected); }
    }
}
void HParticleTrackSorter::cleanUp(Bool_t final)
{
    // delete temporary candidate objects from
    // heap and clear the candidate vector. clear
    // all index vectors.

    clearVector(all_candidates);

    if(kIgnorePreviousIndex || final){
	index_RICH    .clear();
	index_InnerMDC.clear();
	index_OuterMDC.clear();
	index_SHOWER  .clear();
	index_TOFHit  .clear();
	index_TOFClst .clear();
	index_RPC     .clear();
    }
}

Int_t  HParticleTrackSorter::fillAndSetFlags()
{
    // clean up the memory and index vectors and fill the
    // candidate vector new. Flag the kIsAcceptedXXX criteria
    // and kIsDoubleHitXXX criteria of the HParticleCand objects of
    // the selection (skipp objects which are marked kIsUsed and
    // kIsRejected) returns the number of candidates.


    if(!pParticleCandCat) return 0; // nothing to do
    fill_Iteration++;

    cleanUp(kFALSE);
    //------------------------------------------------------------------
    // fill HParticleCand Info to temorary objects
    // which will be jused for sorting
    // all flags in HParticleCand will be reset (besides kIsRejected)
    fillInput(all_candidates);


    // list of thing to to
    // for each Double_t hit
    Particle::ESwitch myIndex  [] = {
	Particle::kIsIndexRICHSorter,
	Particle::kIsIndexInnerMDCSorter,
	Particle::kIsIndexOuterMDCSorter,
	Particle::kIsIndexTOFSorter,
	Particle::kIsIndexSHOWERSorter,
        Particle::kIsIndexRPCSorter
    };

    // and for different quality criteria
    Particle::ESwitch myQuality[] = {
        Particle::kIsBestHitRICHSorter,
        Particle::kIsBestHitInnerMDCSorter,
        Particle::kIsBestHitOuterMDCSorter,
        Particle::kIsBestHitMETASorter,
        Particle::kIsBestRKSorter,
	Particle::kIsBestRKRKMETASorter,
	Particle::kIsBestRKRKMETARadiusSorter,
        Particle::kIsBestUserSorter
   };

    kSwitchIndex   = Particle::kIsIndexRICHSorter;
    kSwitchQuality = Particle::kIsBestHitRICHSorter;

    // check all objects for validity with respect to
    // quality criteria and flag them in HParticleCand
    for(UInt_t i = 0; i < (sizeof(myQuality) / sizeof(Int_t)); i ++){
	flagAccepted(all_candidates,myQuality[i]);
    }

    // vector is resorted by index! flag all Double_t used
    // detector hits in HParticleCand
    for(UInt_t i = 0; i < (sizeof(myIndex) / sizeof(Int_t)); i ++){
        flagDouble(all_candidates,myIndex[i]);
    }

    return all_candidates.empty() ? 0 : all_candidates.size();
}

void  HParticleTrackSorter::backupFlags(Bool_t onlyFlags)
{
    // backup the flags of HParticleCand objects. The flags cand be
    // restored to the objects by calling restoreFlags. The flags
    // which have been stored previously will be cleared.

    if(!onlyFlags){
	setup_backup.copyFromStatic();
    }
    old_flags.clear();
    old_flags_emc.clear();

    if(!iterParticleCandCat) return;

    HParticleCand* pCand ;
    iterParticleCandCat->Reset();

    while ((pCand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0 )
    {
         old_flags.push_back(pCand->getFlagField());
    }
    if(pEmcClusterCat){
	HEmcCluster* emc = 0;
	Int_t n = pEmcClusterCat->getEntries();
	for(Int_t i = 0; i < n; i ++){
	    emc = (HEmcCluster*) pEmcClusterCat->getObject(i);
	    if(emc) {
		Int_t fl = emc -> isUsedInParticleCand();
		old_flags_emc.push_back(fl);
	    }
	}
    }

}

Bool_t  HParticleTrackSorter::restoreFlags(Bool_t onlyFlags)
{
    // restore the flags of HParticleCand objects since the last call of
    // backupFlags(). Returns kFALSE if size of vector old_flags does
    // not match the number of objects in HParticleCand category

    if(!onlyFlags){
	setup_backup.copyToStatic();
    }

    if(!iterParticleCandCat) return kFALSE;

    HParticleCand* pCand ;
    iterParticleCandCat->Reset();

    vector<Int_t>::iterator iter;
    iter = old_flags.begin();
    while ((pCand = dynamic_cast<HParticleCand*>(iterParticleCandCat->Next())) != 0 )
    {
	if(iter != old_flags.end()){
	    pCand->setFlagField(*iter);
	    iter ++;
	} else {
            Error("restoreFlags()","Size of backuped flags does not match number of HParticleCand objects!");
	    return kFALSE;
	}
    }

    if(pEmcClusterCat){
	HEmcCluster* emc = 0;
	Int_t n = pEmcClusterCat->getEntries();
	for(Int_t i = 0; i < n; i ++){
	    emc = (HEmcCluster*) pEmcClusterCat->getObject(i);
	    if(emc) {
		emc -> unsetIsUsedInParticleCand();
                if(old_flags_emc[i] == 1 ) emc -> isUsedInParticleCand();
	    }
	}
    }

    return kTRUE;
}
void  HParticleTrackSorter::backupSetup()
{
    // backup the setup of HParticleTrackSorter. The setup cand be
    // restored to the objects by calling restoreSetup.

    setup_backup.copyFromStatic();
}

Bool_t  HParticleTrackSorter::restoreSetup()
{
    // restore the setup of HParticleTrackSorter since the last call of
    // backupSetup().

    setup_backup.copyToStatic();

    return kTRUE;
}

Int_t HParticleTrackSorter::fill(Bool_t (*function)(HParticleCand* ))
{
    // fill the input candidate vector. The kIsRejected flag and
    // lower flags (below kIsLepton) is reseted for all objects
    // not marked as kIsUsed. The filling is done twice once for
    // all objects and second after the application of the user
    // provide test function. Returns the number of candidates after
    // the selection. If kDebug is set and printLevel >= 1 the number
    // of candidates before and after selection will be printed.

    resetFlags(kTRUE,kTRUE,kFALSE,kFALSE);  // need to reset reject/double hit bits (already used objects will stay)!
    Int_t first  = fillAndSetFlags();       // fill all candidates
    selection(function);                    // flag candidates rejected if the they fail the test
    Int_t second = fillAndSetFlags();       // fill only good candidates
    if(kDebug && printLevel >= 1) {
	cout<<"fill(): "
	    <<"s: "                <<selectBest_Iteration
	    <<", f: "              <<fill_Iteration
	    <<", kSwitchParticle: "<<kSwitchParticle
	    <<", nCands first: "   <<first
	    <<", nCands second :"  <<second
	    <<endl;
    }
    return all_candidates.empty() ?  0 : all_candidates.size();
}

Int_t HParticleTrackSorter::selectBest(Particle::ESwitch byQuality, Int_t byParticle)
{
    // perfomes the selection of the best tracks (with respect to the chosen quality
    // criteria), which makes shure, that each detector hit is used only once. If
    // In case Int_t byParticle == kIsLepton the RICH index
    // is taken into account, otherwise not. For leptons the kIsLepton bit is flagged
    // in addition. Individual detector hits can be ignored in the selection procedure with
    // the setIgnoreXXX(kTRUE) functions (by default all hits are taken into account).
    // Returns the number of selected best HParticleCand objects.
    // if kDebug is set and  printLevel >= 1 the number of sected objects and
    // the statistics of the rejection on Double_t hits i the snhgle detectors is printed.
    // if printLevel >=3 the properties of the candidates (accepted and rejected) are
    // printed.

    if(fEvent->getHeader()->getEventSeqNumber() != currentEvent){
	currentEvent         = fEvent->getHeader()->getEventSeqNumber();
        selectBest_Iteration = 0;
    }
    selectBest_Iteration++;

    kSwitchQuality  = byQuality;
    kSwitchParticle = byParticle;
    Int_t ct        = 0;
    Int_t ctdouble[4] = {0};

    if(kSwitchQuality==Particle::kIsBestUserSorter && pUserSort==NULL){
	::Error ("selectBest()","qualtiy is kIsBestUser, but pUserSort==NULL. No sort done.");
    } else {
	sort(all_candidates.begin(),all_candidates.end(),this->cmpQuality);
    }
    for(UInt_t i = 0; i < all_candidates.size(); ++ i)
    {
        candidateSort* cand = all_candidates[i];

        //-------------------------------------------------------
	// skip the candidates which contain a single detector hit
        // which has been used before
        Int_t isDoubleHit = 0;
        /*
	// disabled :: HParticleCand objects are not created for all comb.
	// of innerMDC  and RICH ! This caused problems with close pairs
	// (2 good tracks matched to same RICH ring).
	if( (kSwitchParticle == Particle::kIsLeptonSorter && !kIgnoreRICH ) ) {
            // look to RICH hits only in case of leptons
	    if(cand->ind_RICH     >= 0 &&
	       (find(index_RICH.begin(),index_RICH.end()        ,cand->ind_RICH)     != index_RICH.end()) ){
		isDoubleHit |= 0x1 << 0;
		ctdouble[0] ++;
	    }
	} */
	if( !kIgnoreInnerMDC ) {
	    if(cand->ind_innerMDC >= 0 &&
	       (find(index_InnerMDC.begin(),index_InnerMDC.end(),cand->ind_innerMDC) != index_InnerMDC.end()) ){
		isDoubleHit |= 0x1 << 1;
		ctdouble[1] ++;
	    }
	}
	if( !kIgnoreOuterMDC ) {
	    if(cand->ind_outerMDC >= 0 &&
	       (find(index_OuterMDC.begin(),index_OuterMDC.end(),cand->ind_outerMDC) != index_OuterMDC.end()) ){
		isDoubleHit |= 0x1 << 2;
		ctdouble[2] ++;
	    }
	}
	if( !kIgnoreMETA ) {
	    if( (cand->selectedMeta == Particle::kTofHit1 || cand->selectedMeta == Particle::kTofHit2)  &&
	       cand->ind_TOF      >= 0 &&
	       (find(index_TOFHit.begin(),index_TOFHit.end()    ,cand->ind_TOF)      != index_TOFHit.end()) ){
		isDoubleHit |= 0x1 << 3;
		ctdouble[3] ++;

	    }

	    if(cand->selectedMeta == Particle::kTofClst && cand->ind_TOF      >= 0 &&
	       (find(index_TOFClst.begin(),index_TOFClst.end()   ,cand->ind_TOF)     != index_TOFClst.end()) ){
		isDoubleHit |= 0x1 << 3;
		ctdouble[3] ++;

	    }
            /*
	    // disabled :: HParticleCand objects are not created for all comb.
	    // of outerMDC  and SHOWER!
	    if(cand->selectedMeta == Particle::kShowerHit && cand->ind_SHOWER   >= 0 &&
	       (find(index_SHOWER.begin(),index_SHOWER.end()    ,cand->ind_SHOWER)   != index_SHOWER.end()) ){
		isDoubleHit |= 0x1 << 3;
		ctdouble[3] ++;
	    }
            */
	    if(cand->selectedMeta == Particle::kRpcClst   && cand->ind_RPC      >= 0 &&
	       (find(index_RPC.begin(),index_RPC.end()          ,cand->ind_RPC)      != index_RPC.end()) ){
		isDoubleHit |= 0x1 << 3;
		ctdouble[3] ++;
	    }
	}
	if(isDoubleHit){
            if(kDebug && printLevel >= 1){
		TString out="";
		for(Int_t j = 3; j >= 0; j --){
		    out += ( (isDoubleHit>>j) & 0x1 );
		}
		cout<<"--selectBest(): skip Double_t hit with pattern : "<<out.Data()<<endl;
		printCand(cand,i);
	    }
	    continue;
	}
        //-------------------------------------------------------

        //-------------------------------------------------------
	// mark the candidate as best and add the single detector hit indices
	// to the lookup

	/*
	// disabled :: HParticleCand objects are not created for all comb.
	// of innerMDC and RICH ! This caused problems with close pairs
	// (2 good tracks matched to same RICH ring).
	if(cand->ind_RICH     >= 0 ){
	    index_RICH.push_back(cand->ind_RICH);
	}
	*/

	if(cand->ind_innerMDC >= 0 ){
	    index_InnerMDC.push_back(cand->ind_innerMDC);
	}
	if(cand->ind_outerMDC >= 0 ){
	    index_OuterMDC.push_back(cand->ind_outerMDC);
	}
	if( (cand->selectedMeta == Particle::kTofHit1 || cand->selectedMeta == Particle::kTofHit2) &&
	   cand->ind_TOF      >= 0 ){
	    index_TOFHit.push_back(cand->ind_TOF);
	}
	if(cand->selectedMeta == Particle::kTofClst &&
	   cand->ind_TOF      >= 0 ){
	    index_TOFClst.push_back(cand->ind_TOF);
	}
        /*
        // disabled :: HParticleCand objects are not created for all comb.
	// of outerMDC and SHOWER !
	if(cand->selectedMeta == Particle::kShowerHit &&
	   cand->ind_SHOWER   >= 0 ){
	    index_SHOWER.push_back(cand->ind_SHOWER);
	}
	*/
	if(cand->selectedMeta == Particle::kRpcClst &&
	   cand->ind_RPC      >= 0 ){
	    index_RPC.push_back(cand->ind_RPC);
	}

	HParticleCand* pcand = dynamic_cast<HParticleCand*>(pParticleCandCat->getObject(cand->ind_Cand));
	if(!pcand) {
	    Error("HParticleTrackSorter::selectBest()","NULL pointer retrieved for HParticleCand at index %i",cand->ind_Cand);
	    continue;
	}
	switch (byQuality)
	{
	case Particle::kIsBestHitRICHSorter:
	    pcand -> setFlagBit(Particle::kIsBestHitRICH);
	    break;
	case Particle::kIsBestHitInnerMDCSorter:
	    pcand -> setFlagBit(Particle::kIsBestHitInnerMDC);
	    break;
	case Particle::kIsBestHitOuterMDCSorter:
	    pcand -> setFlagBit(Particle::kIsBestHitOuterMDC);
	    break;
	case Particle::kIsBestHitMETASorter:
	    pcand -> setFlagBit(Particle::kIsBestHitMETA);
	    break;
	case Particle::kIsBestRKSorter:
	    pcand -> setFlagBit(Particle::kIsBestRK);
	    break;
	case Particle::kIsBestRKRKMETASorter:
	    pcand -> setFlagBit(Particle::kIsBestRKRKMETA);
	    break;
	case Particle::kIsBestRKRKMETARadiusSorter:
	    pcand -> setFlagBit(Particle::kIsBestRKRKMETARadius);
	    break;
	case Particle::kIsBestUserSorter:
	    pcand -> setFlagBit(Particle::kIsBestUser);
	    break;
	default :

	    Error("HParticleTrackSorter::selectBest()","Unknown switch statement : %i!",byQuality);
	    exit(1);
	    break;

	}
        ++ ct;
	pcand ->setFlagBit(Particle::kIsUsed);
	if(byParticle == Particle::kIsLeptonSorter){
	    pcand ->setFlagBit(Particle::kIsLepton);
	}
	if(kDebug && printLevel >= 1){
	    cout<<"++selectBest(): select hit : "<<endl;
	    printCand(cand,i);
	}
	//-------------------------------------------------------

    }
    if(kDebug && printLevel >= 1)
    {
	cout<<"selectBest(): "
	    <<"s: "                  <<selectBest_Iteration
	    <<", f: "                <<fill_Iteration
	    <<", kSwitchParticle: "  <<kSwitchParticle
	    <<", nBest cand: "       <<ct
	    <<", Double_t RICH: "    <<ctdouble[0]
            <<", InnerMDC: "         <<ctdouble[1]
            <<", Double_t OuterMDC: "<<ctdouble[2]
            <<", Double_t META: "    <<ctdouble[3]
	    <<endl;
    }
    fill_Iteration = 0;
    flagEmcClusters();
  return ct;
}

void HParticleTrackSorter::setRICHMatching(Particle::ERichMatch match, Float_t window)
{
    // switch between the matching cut RICH-MDC inside the selectLeptons() function. The method
    // is selected by the options defined in Particle::ERichMatch.
	kSwitchRICHMatching = match;
	fRICHMDCWindow      = window;
	setup.kSwitchRICHMatching = match;
	setup.fRICHMDCWindow      = window;
}


Bool_t HParticleTrackSorter::selectLeptons(HParticleCand* pcand)
{
    // build in selection function for lepton candidates.
    // Requires besides an RICH hit, RK + META and fitted
    // inner+outer segment. The RICH-MDC matching cut depends on the
    // settings of the matching method and window (set by
    // HParticleTrackSorter::setRICHMatching . default is
    // kUseRKRICHWindow, window = 4. degree)


    if(kUseFakeRejection&&pcand->isFakeRejected()) return kFALSE;
    Bool_t test = kFALSE;
    if( pcand->isFlagAND(5,
			 Particle::kIsAcceptedHitRICH,
			 Particle::kIsAcceptedHitInnerMDC,
			 Particle::kIsAcceptedHitOuterMDC,
			 Particle::kIsAcceptedHitMETA,
			 Particle::kIsAcceptedRK
			)
       &&
       pcand->getInnerSegmentChi2() > 0
       &&
       pcand->getChi2()             < 1000      // RK
      ) test = kTRUE;

    if(!test) return kFALSE;

    if(kUseMETAQA && (pcand->getMetaMatchQuality() > fMETAQACut || pcand->getMetaMatchQuality() < 0)) return kFALSE;
    if(kUseYMatching && !HParticleTool::isGoodMetaCell(pcand,fMetaBoundary,kUseYMatchingScaling)) return kFALSE;
    if(kUseBeta && pcand->getBeta()<fBetaLepCut) return kFALSE;
    if(!kIsFieldOn && !isDeltaSegOk(pcand)) return kFALSE;

    if(test)
    {   // rich match only!
	if( kSwitchRICHMatching == Particle::kUseRICHIndexSorter)
	{
	    return kTRUE;
	}
	else if ( kSwitchRICHMatching == Particle::kUseRKRICHCorrelationSorter	&&  pcand -> getRingCorr() )
	{
	    return kTRUE;
	}
	else if ( kSwitchRICHMatching == Particle::kUseRKRICHWindowSorter
		 &&
		 fabs(pcand->getDeltaTheta()) < fRICHMDCWindow
		 &&
		 fabs(pcand->getDeltaPhi())   < fRICHMDCWindow
		)
	{
	    return kTRUE;
	}
	else return kFALSE;

    }
    else return kFALSE;
}
Bool_t HParticleTrackSorter::selectLeptonsNoOuterFit(HParticleCand* pcand)
{
    // build in selection function for lepton candidates.
    // Requires besides an RICH hit, RK + META and fitted
    // inner. Fit for outer segment not reuired. The
    // RICH-MDC matching cut depends on the
    // settings of the matching method and window (set by
    // HParticleTrackSorter::setRICHMatching . default is
    // kUseRKRICHWindow, window = 4. degree)

    if(kUseFakeRejection&&pcand->isFakeRejected()) return kFALSE;
    Bool_t test = kFALSE;
    if( pcand->isFlagAND(4,
			 Particle::kIsAcceptedHitRICH,
			 Particle::kIsAcceptedHitInnerMDC,
			 Particle::kIsAcceptedHitMETA,
			 Particle::kIsAcceptedRK
			)
       &&
       pcand->getInnerSegmentChi2() > 0
       &&
       pcand->getOuterSegInd()      > 0
       &&
       pcand->getChi2()             < 1000      // RK
      ) test = kTRUE;

    if(!test) return kFALSE;

    if(kUseMETAQA && (pcand->getMetaMatchQuality() > fMETAQACut || pcand->getMetaMatchQuality() < 0)) return kFALSE;
    if(kUseYMatching && !HParticleTool::isGoodMetaCell(pcand,fMetaBoundary,kUseYMatchingScaling)) return kFALSE;
    if(kUseBeta && pcand->getBeta()<fBetaLepCut) return kFALSE;
    if(!kIsFieldOn && !isDeltaSegOk(pcand)) return kFALSE;

    if(test)
    {   // rich match only!
	if( kSwitchRICHMatching == Particle::kUseRICHIndexSorter)
	{
	    return kTRUE;
	}
	else if ( kSwitchRICHMatching == Particle::kUseRKRICHCorrelationSorter	&&  pcand -> getRingCorr() )
	{
	    return kTRUE;
	}
	else if ( kSwitchRICHMatching == Particle::kUseRKRICHWindowSorter
		 &&
		 fabs(pcand->getDeltaTheta()) < fRICHMDCWindow
		 &&
		 fabs(pcand->getDeltaPhi())   < fRICHMDCWindow
		)
	{
	    return kTRUE;
	}
	else return kFALSE;
    }
    else return kFALSE;
}
Bool_t HParticleTrackSorter::selectHadrons(HParticleCand* pcand )
{
    // build in selection function for hadron candidates.
    // Requires besides RK + META and fitted inner+outer MDC.

    if(kUseFakeRejection&&pcand->isFakeRejected()) return kFALSE;
    Bool_t test=kFALSE;
    if( pcand->isFlagAND(4,
			 Particle::kIsAcceptedHitInnerMDC,
			 Particle::kIsAcceptedHitOuterMDC,
			 Particle::kIsAcceptedHitMETA,
			 Particle::kIsAcceptedRK
			)
       &&
       pcand->getInnerSegmentChi2() > 0
       &&
       pcand->getChi2()             < 1000
      ) test=kTRUE;
    if(!test) return kFALSE;

    if(kUseMETAQA && (pcand->getMetaMatchQuality() > fMETAQACut || pcand->getMetaMatchQuality() < 0)) return kFALSE;
    if(kUseYMatching && !HParticleTool::isGoodMetaCell(pcand,fMetaBoundary,kUseYMatchingScaling)) return kFALSE;
    if(kUseBeta && (pcand->getBeta()<0 || pcand->getTof()>fTOFMAXCut)) return kFALSE;
    if(!kIsFieldOn && !isDeltaSegOk(pcand)) return kFALSE;

    return test;
}
Bool_t HParticleTrackSorter::selectHadronsNoOuterFit(HParticleCand* pcand )
{
    // build in selection function for hadron candidates.
    // Requires besides RK + META and fitted inner MDC. Fit
    // for outer MDC not required.

    if(kUseFakeRejection&&pcand->isFakeRejected()) return kFALSE;
    Bool_t test=kFALSE;
    if( pcand->isFlagAND(3,
			 Particle::kIsAcceptedHitInnerMDC,
			 Particle::kIsAcceptedHitMETA,
			 Particle::kIsAcceptedRK
			)
       &&
       pcand->getInnerSegmentChi2() > 0
       &&
       pcand->getOuterSegInd()      > 0
       &&
       pcand->getChi2()             < 1000
      ) test=kTRUE;

    if(!test) return kFALSE;

    if(kUseMETAQA && (pcand->getMetaMatchQuality() > fMETAQACut || pcand->getMetaMatchQuality() < 0)) return kFALSE;
    if(kUseYMatching && !HParticleTool::isGoodMetaCell(pcand,fMetaBoundary,kUseYMatchingScaling)) return kFALSE;
    if(kUseBeta && (pcand->getBeta()<0 || pcand->getTof()>fTOFMAXCut)) return kFALSE;
    if(!kIsFieldOn && !isDeltaSegOk(pcand)) return kFALSE;

    return test;
}



void sorter_setup::copyFromStatic()
{

	kSwitchRICHMatching  = HParticleTrackSorter::kSwitchRICHMatching;
	fRICHMDCWindow       = HParticleTrackSorter::fRICHMDCWindow;
	kIsFieldOn           = HParticleTrackSorter::kIsFieldOn;
	kIgnoreRICH          = HParticleTrackSorter::kIgnoreRICH;
	kIgnoreInnerMDC      = HParticleTrackSorter::kIgnoreInnerMDC;
	kIgnoreOuterMDC      = HParticleTrackSorter::kIgnoreOuterMDC;
	kIgnoreMETA          = HParticleTrackSorter::kIgnoreMETA;
	kIgnorePreviousIndex = HParticleTrackSorter::kIgnorePreviousIndex;
	pUserSort            = HParticleTrackSorter::pUserSort;
	kUseYMatching        = HParticleTrackSorter::kUseYMatching;
	kUseYMatchingScaling = HParticleTrackSorter::kUseYMatchingScaling;
	fMetaBoundary        = HParticleTrackSorter::fMetaBoundary;
	kUseBeta             = HParticleTrackSorter::kUseBeta;
	fBetaLepCut          = HParticleTrackSorter::fBetaLepCut;
	kUseFakeRejection    = HParticleTrackSorter::kUseFakeRejection;
	kUseMETAQA           = HParticleTrackSorter::kUseMETAQA;
	fMETAQACut           = HParticleTrackSorter::fMETAQACut;
	fTOFMAXCut           = HParticleTrackSorter::fTOFMAXCut;
	fDeltaSegCut         = HParticleTrackSorter::fDeltaSegCut;
};

void sorter_setup::copyToStatic()
{

	HParticleTrackSorter::kSwitchRICHMatching  = kSwitchRICHMatching ;
	HParticleTrackSorter::fRICHMDCWindow       = fRICHMDCWindow;
	HParticleTrackSorter::kIsFieldOn           = kIsFieldOn;
	HParticleTrackSorter::kIgnoreRICH          = kIgnoreRICH;
	HParticleTrackSorter::kIgnoreInnerMDC      = kIgnoreInnerMDC;
	HParticleTrackSorter::kIgnoreOuterMDC      = kIgnoreOuterMDC;
	HParticleTrackSorter::kIgnoreMETA          = kIgnoreMETA;
	HParticleTrackSorter::kIgnorePreviousIndex = kIgnorePreviousIndex;
	HParticleTrackSorter::pUserSort            = pUserSort;
	HParticleTrackSorter::kUseYMatching        = kUseYMatching;
	HParticleTrackSorter::kUseYMatchingScaling = kUseYMatchingScaling;
	HParticleTrackSorter::fMetaBoundary        = fMetaBoundary;
	HParticleTrackSorter::kUseBeta             = kUseBeta;
	HParticleTrackSorter::fBetaLepCut          = fBetaLepCut;
	HParticleTrackSorter::kUseFakeRejection    = kUseFakeRejection;
	HParticleTrackSorter::kUseMETAQA           = kUseMETAQA;
	HParticleTrackSorter::fMETAQACut           = fMETAQACut;
	HParticleTrackSorter::fTOFMAXCut           = fTOFMAXCut;
	HParticleTrackSorter::fDeltaSegCut         = fDeltaSegCut;
};


 hparticletracksorter.cc:1
 hparticletracksorter.cc:2
 hparticletracksorter.cc:3
 hparticletracksorter.cc:4
 hparticletracksorter.cc:5
 hparticletracksorter.cc:6
 hparticletracksorter.cc:7
 hparticletracksorter.cc:8
 hparticletracksorter.cc:9
 hparticletracksorter.cc:10
 hparticletracksorter.cc:11
 hparticletracksorter.cc:12
 hparticletracksorter.cc:13
 hparticletracksorter.cc:14
 hparticletracksorter.cc:15
 hparticletracksorter.cc:16
 hparticletracksorter.cc:17
 hparticletracksorter.cc:18
 hparticletracksorter.cc:19
 hparticletracksorter.cc:20
 hparticletracksorter.cc:21
 hparticletracksorter.cc:22
 hparticletracksorter.cc:23
 hparticletracksorter.cc:24
 hparticletracksorter.cc:25
 hparticletracksorter.cc:26
 hparticletracksorter.cc:27
 hparticletracksorter.cc:28
 hparticletracksorter.cc:29
 hparticletracksorter.cc:30
 hparticletracksorter.cc:31
 hparticletracksorter.cc:32
 hparticletracksorter.cc:33
 hparticletracksorter.cc:34
 hparticletracksorter.cc:35
 hparticletracksorter.cc:36
 hparticletracksorter.cc:37
 hparticletracksorter.cc:38
 hparticletracksorter.cc:39
 hparticletracksorter.cc:40
 hparticletracksorter.cc:41
 hparticletracksorter.cc:42
 hparticletracksorter.cc:43
 hparticletracksorter.cc:44
 hparticletracksorter.cc:45
 hparticletracksorter.cc:46
 hparticletracksorter.cc:47
 hparticletracksorter.cc:48
 hparticletracksorter.cc:49
 hparticletracksorter.cc:50
 hparticletracksorter.cc:51
 hparticletracksorter.cc:52
 hparticletracksorter.cc:53
 hparticletracksorter.cc:54
 hparticletracksorter.cc:55
 hparticletracksorter.cc:56
 hparticletracksorter.cc:57
 hparticletracksorter.cc:58
 hparticletracksorter.cc:59
 hparticletracksorter.cc:60
 hparticletracksorter.cc:61
 hparticletracksorter.cc:62
 hparticletracksorter.cc:63
 hparticletracksorter.cc:64
 hparticletracksorter.cc:65
 hparticletracksorter.cc:66
 hparticletracksorter.cc:67
 hparticletracksorter.cc:68
 hparticletracksorter.cc:69
 hparticletracksorter.cc:70
 hparticletracksorter.cc:71
 hparticletracksorter.cc:72
 hparticletracksorter.cc:73
 hparticletracksorter.cc:74
 hparticletracksorter.cc:75
 hparticletracksorter.cc:76
 hparticletracksorter.cc:77
 hparticletracksorter.cc:78
 hparticletracksorter.cc:79
 hparticletracksorter.cc:80
 hparticletracksorter.cc:81
 hparticletracksorter.cc:82
 hparticletracksorter.cc:83
 hparticletracksorter.cc:84
 hparticletracksorter.cc:85
 hparticletracksorter.cc:86
 hparticletracksorter.cc:87
 hparticletracksorter.cc:88
 hparticletracksorter.cc:89
 hparticletracksorter.cc:90
 hparticletracksorter.cc:91
 hparticletracksorter.cc:92
 hparticletracksorter.cc:93
 hparticletracksorter.cc:94
 hparticletracksorter.cc:95
 hparticletracksorter.cc:96
 hparticletracksorter.cc:97
 hparticletracksorter.cc:98
 hparticletracksorter.cc:99
 hparticletracksorter.cc:100
 hparticletracksorter.cc:101
 hparticletracksorter.cc:102
 hparticletracksorter.cc:103
 hparticletracksorter.cc:104
 hparticletracksorter.cc:105
 hparticletracksorter.cc:106
 hparticletracksorter.cc:107
 hparticletracksorter.cc:108
 hparticletracksorter.cc:109
 hparticletracksorter.cc:110
 hparticletracksorter.cc:111
 hparticletracksorter.cc:112
 hparticletracksorter.cc:113
 hparticletracksorter.cc:114
 hparticletracksorter.cc:115
 hparticletracksorter.cc:116
 hparticletracksorter.cc:117
 hparticletracksorter.cc:118
 hparticletracksorter.cc:119
 hparticletracksorter.cc:120
 hparticletracksorter.cc:121
 hparticletracksorter.cc:122
 hparticletracksorter.cc:123
 hparticletracksorter.cc:124
 hparticletracksorter.cc:125
 hparticletracksorter.cc:126
 hparticletracksorter.cc:127
 hparticletracksorter.cc:128
 hparticletracksorter.cc:129
 hparticletracksorter.cc:130
 hparticletracksorter.cc:131
 hparticletracksorter.cc:132
 hparticletracksorter.cc:133
 hparticletracksorter.cc:134
 hparticletracksorter.cc:135
 hparticletracksorter.cc:136
 hparticletracksorter.cc:137
 hparticletracksorter.cc:138
 hparticletracksorter.cc:139
 hparticletracksorter.cc:140
 hparticletracksorter.cc:141
 hparticletracksorter.cc:142
 hparticletracksorter.cc:143
 hparticletracksorter.cc:144
 hparticletracksorter.cc:145
 hparticletracksorter.cc:146
 hparticletracksorter.cc:147
 hparticletracksorter.cc:148
 hparticletracksorter.cc:149
 hparticletracksorter.cc:150
 hparticletracksorter.cc:151
 hparticletracksorter.cc:152
 hparticletracksorter.cc:153
 hparticletracksorter.cc:154
 hparticletracksorter.cc:155
 hparticletracksorter.cc:156
 hparticletracksorter.cc:157
 hparticletracksorter.cc:158
 hparticletracksorter.cc:159
 hparticletracksorter.cc:160
 hparticletracksorter.cc:161
 hparticletracksorter.cc:162
 hparticletracksorter.cc:163
 hparticletracksorter.cc:164
 hparticletracksorter.cc:165
 hparticletracksorter.cc:166
 hparticletracksorter.cc:167
 hparticletracksorter.cc:168
 hparticletracksorter.cc:169
 hparticletracksorter.cc:170
 hparticletracksorter.cc:171
 hparticletracksorter.cc:172
 hparticletracksorter.cc:173
 hparticletracksorter.cc:174
 hparticletracksorter.cc:175
 hparticletracksorter.cc:176
 hparticletracksorter.cc:177
 hparticletracksorter.cc:178
 hparticletracksorter.cc:179
 hparticletracksorter.cc:180
 hparticletracksorter.cc:181
 hparticletracksorter.cc:182
 hparticletracksorter.cc:183
 hparticletracksorter.cc:184
 hparticletracksorter.cc:185
 hparticletracksorter.cc:186
 hparticletracksorter.cc:187
 hparticletracksorter.cc:188
 hparticletracksorter.cc:189
 hparticletracksorter.cc:190
 hparticletracksorter.cc:191
 hparticletracksorter.cc:192
 hparticletracksorter.cc:193
 hparticletracksorter.cc:194
 hparticletracksorter.cc:195
 hparticletracksorter.cc:196
 hparticletracksorter.cc:197
 hparticletracksorter.cc:198
 hparticletracksorter.cc:199
 hparticletracksorter.cc:200
 hparticletracksorter.cc:201
 hparticletracksorter.cc:202
 hparticletracksorter.cc:203
 hparticletracksorter.cc:204
 hparticletracksorter.cc:205
 hparticletracksorter.cc:206
 hparticletracksorter.cc:207
 hparticletracksorter.cc:208
 hparticletracksorter.cc:209
 hparticletracksorter.cc:210
 hparticletracksorter.cc:211
 hparticletracksorter.cc:212
 hparticletracksorter.cc:213
 hparticletracksorter.cc:214
 hparticletracksorter.cc:215
 hparticletracksorter.cc:216
 hparticletracksorter.cc:217
 hparticletracksorter.cc:218
 hparticletracksorter.cc:219
 hparticletracksorter.cc:220
 hparticletracksorter.cc:221
 hparticletracksorter.cc:222
 hparticletracksorter.cc:223
 hparticletracksorter.cc:224
 hparticletracksorter.cc:225
 hparticletracksorter.cc:226
 hparticletracksorter.cc:227
 hparticletracksorter.cc:228
 hparticletracksorter.cc:229
 hparticletracksorter.cc:230
 hparticletracksorter.cc:231
 hparticletracksorter.cc:232
 hparticletracksorter.cc:233
 hparticletracksorter.cc:234
 hparticletracksorter.cc:235
 hparticletracksorter.cc:236
 hparticletracksorter.cc:237
 hparticletracksorter.cc:238
 hparticletracksorter.cc:239
 hparticletracksorter.cc:240
 hparticletracksorter.cc:241
 hparticletracksorter.cc:242
 hparticletracksorter.cc:243
 hparticletracksorter.cc:244
 hparticletracksorter.cc:245
 hparticletracksorter.cc:246
 hparticletracksorter.cc:247
 hparticletracksorter.cc:248
 hparticletracksorter.cc:249
 hparticletracksorter.cc:250
 hparticletracksorter.cc:251
 hparticletracksorter.cc:252
 hparticletracksorter.cc:253
 hparticletracksorter.cc:254
 hparticletracksorter.cc:255
 hparticletracksorter.cc:256
 hparticletracksorter.cc:257
 hparticletracksorter.cc:258
 hparticletracksorter.cc:259
 hparticletracksorter.cc:260
 hparticletracksorter.cc:261
 hparticletracksorter.cc:262
 hparticletracksorter.cc:263
 hparticletracksorter.cc:264
 hparticletracksorter.cc:265
 hparticletracksorter.cc:266
 hparticletracksorter.cc:267
 hparticletracksorter.cc:268
 hparticletracksorter.cc:269
 hparticletracksorter.cc:270
 hparticletracksorter.cc:271
 hparticletracksorter.cc:272
 hparticletracksorter.cc:273
 hparticletracksorter.cc:274
 hparticletracksorter.cc:275
 hparticletracksorter.cc:276
 hparticletracksorter.cc:277
 hparticletracksorter.cc:278
 hparticletracksorter.cc:279
 hparticletracksorter.cc:280
 hparticletracksorter.cc:281
 hparticletracksorter.cc:282
 hparticletracksorter.cc:283
 hparticletracksorter.cc:284
 hparticletracksorter.cc:285
 hparticletracksorter.cc:286
 hparticletracksorter.cc:287
 hparticletracksorter.cc:288
 hparticletracksorter.cc:289
 hparticletracksorter.cc:290
 hparticletracksorter.cc:291
 hparticletracksorter.cc:292
 hparticletracksorter.cc:293
 hparticletracksorter.cc:294
 hparticletracksorter.cc:295
 hparticletracksorter.cc:296
 hparticletracksorter.cc:297
 hparticletracksorter.cc:298
 hparticletracksorter.cc:299
 hparticletracksorter.cc:300
 hparticletracksorter.cc:301
 hparticletracksorter.cc:302
 hparticletracksorter.cc:303
 hparticletracksorter.cc:304
 hparticletracksorter.cc:305
 hparticletracksorter.cc:306
 hparticletracksorter.cc:307
 hparticletracksorter.cc:308
 hparticletracksorter.cc:309
 hparticletracksorter.cc:310
 hparticletracksorter.cc:311
 hparticletracksorter.cc:312
 hparticletracksorter.cc:313
 hparticletracksorter.cc:314
 hparticletracksorter.cc:315
 hparticletracksorter.cc:316
 hparticletracksorter.cc:317
 hparticletracksorter.cc:318
 hparticletracksorter.cc:319
 hparticletracksorter.cc:320
 hparticletracksorter.cc:321
 hparticletracksorter.cc:322
 hparticletracksorter.cc:323
 hparticletracksorter.cc:324
 hparticletracksorter.cc:325
 hparticletracksorter.cc:326
 hparticletracksorter.cc:327
 hparticletracksorter.cc:328
 hparticletracksorter.cc:329
 hparticletracksorter.cc:330
 hparticletracksorter.cc:331
 hparticletracksorter.cc:332
 hparticletracksorter.cc:333
 hparticletracksorter.cc:334
 hparticletracksorter.cc:335
 hparticletracksorter.cc:336
 hparticletracksorter.cc:337
 hparticletracksorter.cc:338
 hparticletracksorter.cc:339
 hparticletracksorter.cc:340
 hparticletracksorter.cc:341
 hparticletracksorter.cc:342
 hparticletracksorter.cc:343
 hparticletracksorter.cc:344
 hparticletracksorter.cc:345
 hparticletracksorter.cc:346
 hparticletracksorter.cc:347
 hparticletracksorter.cc:348
 hparticletracksorter.cc:349
 hparticletracksorter.cc:350
 hparticletracksorter.cc:351
 hparticletracksorter.cc:352
 hparticletracksorter.cc:353
 hparticletracksorter.cc:354
 hparticletracksorter.cc:355
 hparticletracksorter.cc:356
 hparticletracksorter.cc:357
 hparticletracksorter.cc:358
 hparticletracksorter.cc:359
 hparticletracksorter.cc:360
 hparticletracksorter.cc:361
 hparticletracksorter.cc:362
 hparticletracksorter.cc:363
 hparticletracksorter.cc:364
 hparticletracksorter.cc:365
 hparticletracksorter.cc:366
 hparticletracksorter.cc:367
 hparticletracksorter.cc:368
 hparticletracksorter.cc:369
 hparticletracksorter.cc:370
 hparticletracksorter.cc:371
 hparticletracksorter.cc:372
 hparticletracksorter.cc:373
 hparticletracksorter.cc:374
 hparticletracksorter.cc:375
 hparticletracksorter.cc:376
 hparticletracksorter.cc:377
 hparticletracksorter.cc:378
 hparticletracksorter.cc:379
 hparticletracksorter.cc:380
 hparticletracksorter.cc:381
 hparticletracksorter.cc:382
 hparticletracksorter.cc:383
 hparticletracksorter.cc:384
 hparticletracksorter.cc:385
 hparticletracksorter.cc:386
 hparticletracksorter.cc:387
 hparticletracksorter.cc:388
 hparticletracksorter.cc:389
 hparticletracksorter.cc:390
 hparticletracksorter.cc:391
 hparticletracksorter.cc:392
 hparticletracksorter.cc:393
 hparticletracksorter.cc:394
 hparticletracksorter.cc:395
 hparticletracksorter.cc:396
 hparticletracksorter.cc:397
 hparticletracksorter.cc:398
 hparticletracksorter.cc:399
 hparticletracksorter.cc:400
 hparticletracksorter.cc:401
 hparticletracksorter.cc:402
 hparticletracksorter.cc:403
 hparticletracksorter.cc:404
 hparticletracksorter.cc:405
 hparticletracksorter.cc:406
 hparticletracksorter.cc:407
 hparticletracksorter.cc:408
 hparticletracksorter.cc:409
 hparticletracksorter.cc:410
 hparticletracksorter.cc:411
 hparticletracksorter.cc:412
 hparticletracksorter.cc:413
 hparticletracksorter.cc:414
 hparticletracksorter.cc:415
 hparticletracksorter.cc:416
 hparticletracksorter.cc:417
 hparticletracksorter.cc:418
 hparticletracksorter.cc:419
 hparticletracksorter.cc:420
 hparticletracksorter.cc:421
 hparticletracksorter.cc:422
 hparticletracksorter.cc:423
 hparticletracksorter.cc:424
 hparticletracksorter.cc:425
 hparticletracksorter.cc:426
 hparticletracksorter.cc:427
 hparticletracksorter.cc:428
 hparticletracksorter.cc:429
 hparticletracksorter.cc:430
 hparticletracksorter.cc:431
 hparticletracksorter.cc:432
 hparticletracksorter.cc:433
 hparticletracksorter.cc:434
 hparticletracksorter.cc:435
 hparticletracksorter.cc:436
 hparticletracksorter.cc:437
 hparticletracksorter.cc:438
 hparticletracksorter.cc:439
 hparticletracksorter.cc:440
 hparticletracksorter.cc:441
 hparticletracksorter.cc:442
 hparticletracksorter.cc:443
 hparticletracksorter.cc:444
 hparticletracksorter.cc:445
 hparticletracksorter.cc:446
 hparticletracksorter.cc:447
 hparticletracksorter.cc:448
 hparticletracksorter.cc:449
 hparticletracksorter.cc:450
 hparticletracksorter.cc:451
 hparticletracksorter.cc:452
 hparticletracksorter.cc:453
 hparticletracksorter.cc:454
 hparticletracksorter.cc:455
 hparticletracksorter.cc:456
 hparticletracksorter.cc:457
 hparticletracksorter.cc:458
 hparticletracksorter.cc:459
 hparticletracksorter.cc:460
 hparticletracksorter.cc:461
 hparticletracksorter.cc:462
 hparticletracksorter.cc:463
 hparticletracksorter.cc:464
 hparticletracksorter.cc:465
 hparticletracksorter.cc:466
 hparticletracksorter.cc:467
 hparticletracksorter.cc:468
 hparticletracksorter.cc:469
 hparticletracksorter.cc:470
 hparticletracksorter.cc:471
 hparticletracksorter.cc:472
 hparticletracksorter.cc:473
 hparticletracksorter.cc:474
 hparticletracksorter.cc:475
 hparticletracksorter.cc:476
 hparticletracksorter.cc:477
 hparticletracksorter.cc:478
 hparticletracksorter.cc:479
 hparticletracksorter.cc:480
 hparticletracksorter.cc:481
 hparticletracksorter.cc:482
 hparticletracksorter.cc:483
 hparticletracksorter.cc:484
 hparticletracksorter.cc:485
 hparticletracksorter.cc:486
 hparticletracksorter.cc:487
 hparticletracksorter.cc:488
 hparticletracksorter.cc:489
 hparticletracksorter.cc:490
 hparticletracksorter.cc:491
 hparticletracksorter.cc:492
 hparticletracksorter.cc:493
 hparticletracksorter.cc:494
 hparticletracksorter.cc:495
 hparticletracksorter.cc:496
 hparticletracksorter.cc:497
 hparticletracksorter.cc:498
 hparticletracksorter.cc:499
 hparticletracksorter.cc:500
 hparticletracksorter.cc:501
 hparticletracksorter.cc:502
 hparticletracksorter.cc:503
 hparticletracksorter.cc:504
 hparticletracksorter.cc:505
 hparticletracksorter.cc:506
 hparticletracksorter.cc:507
 hparticletracksorter.cc:508
 hparticletracksorter.cc:509
 hparticletracksorter.cc:510
 hparticletracksorter.cc:511
 hparticletracksorter.cc:512
 hparticletracksorter.cc:513
 hparticletracksorter.cc:514
 hparticletracksorter.cc:515
 hparticletracksorter.cc:516
 hparticletracksorter.cc:517
 hparticletracksorter.cc:518
 hparticletracksorter.cc:519
 hparticletracksorter.cc:520
 hparticletracksorter.cc:521
 hparticletracksorter.cc:522
 hparticletracksorter.cc:523
 hparticletracksorter.cc:524
 hparticletracksorter.cc:525
 hparticletracksorter.cc:526
 hparticletracksorter.cc:527
 hparticletracksorter.cc:528
 hparticletracksorter.cc:529
 hparticletracksorter.cc:530
 hparticletracksorter.cc:531
 hparticletracksorter.cc:532
 hparticletracksorter.cc:533
 hparticletracksorter.cc:534
 hparticletracksorter.cc:535
 hparticletracksorter.cc:536
 hparticletracksorter.cc:537
 hparticletracksorter.cc:538
 hparticletracksorter.cc:539
 hparticletracksorter.cc:540
 hparticletracksorter.cc:541
 hparticletracksorter.cc:542
 hparticletracksorter.cc:543
 hparticletracksorter.cc:544
 hparticletracksorter.cc:545
 hparticletracksorter.cc:546
 hparticletracksorter.cc:547
 hparticletracksorter.cc:548
 hparticletracksorter.cc:549
 hparticletracksorter.cc:550
 hparticletracksorter.cc:551
 hparticletracksorter.cc:552
 hparticletracksorter.cc:553
 hparticletracksorter.cc:554
 hparticletracksorter.cc:555
 hparticletracksorter.cc:556
 hparticletracksorter.cc:557
 hparticletracksorter.cc:558
 hparticletracksorter.cc:559
 hparticletracksorter.cc:560
 hparticletracksorter.cc:561
 hparticletracksorter.cc:562
 hparticletracksorter.cc:563
 hparticletracksorter.cc:564
 hparticletracksorter.cc:565
 hparticletracksorter.cc:566
 hparticletracksorter.cc:567
 hparticletracksorter.cc:568
 hparticletracksorter.cc:569
 hparticletracksorter.cc:570
 hparticletracksorter.cc:571
 hparticletracksorter.cc:572
 hparticletracksorter.cc:573
 hparticletracksorter.cc:574
 hparticletracksorter.cc:575
 hparticletracksorter.cc:576
 hparticletracksorter.cc:577
 hparticletracksorter.cc:578
 hparticletracksorter.cc:579
 hparticletracksorter.cc:580
 hparticletracksorter.cc:581
 hparticletracksorter.cc:582
 hparticletracksorter.cc:583
 hparticletracksorter.cc:584
 hparticletracksorter.cc:585
 hparticletracksorter.cc:586
 hparticletracksorter.cc:587
 hparticletracksorter.cc:588
 hparticletracksorter.cc:589
 hparticletracksorter.cc:590
 hparticletracksorter.cc:591
 hparticletracksorter.cc:592
 hparticletracksorter.cc:593
 hparticletracksorter.cc:594
 hparticletracksorter.cc:595
 hparticletracksorter.cc:596
 hparticletracksorter.cc:597
 hparticletracksorter.cc:598
 hparticletracksorter.cc:599
 hparticletracksorter.cc:600
 hparticletracksorter.cc:601
 hparticletracksorter.cc:602
 hparticletracksorter.cc:603
 hparticletracksorter.cc:604
 hparticletracksorter.cc:605
 hparticletracksorter.cc:606
 hparticletracksorter.cc:607
 hparticletracksorter.cc:608
 hparticletracksorter.cc:609
 hparticletracksorter.cc:610
 hparticletracksorter.cc:611
 hparticletracksorter.cc:612
 hparticletracksorter.cc:613
 hparticletracksorter.cc:614
 hparticletracksorter.cc:615
 hparticletracksorter.cc:616
 hparticletracksorter.cc:617
 hparticletracksorter.cc:618
 hparticletracksorter.cc:619
 hparticletracksorter.cc:620
 hparticletracksorter.cc:621
 hparticletracksorter.cc:622
 hparticletracksorter.cc:623
 hparticletracksorter.cc:624
 hparticletracksorter.cc:625
 hparticletracksorter.cc:626
 hparticletracksorter.cc:627
 hparticletracksorter.cc:628
 hparticletracksorter.cc:629
 hparticletracksorter.cc:630
 hparticletracksorter.cc:631
 hparticletracksorter.cc:632
 hparticletracksorter.cc:633
 hparticletracksorter.cc:634
 hparticletracksorter.cc:635
 hparticletracksorter.cc:636
 hparticletracksorter.cc:637
 hparticletracksorter.cc:638
 hparticletracksorter.cc:639
 hparticletracksorter.cc:640
 hparticletracksorter.cc:641
 hparticletracksorter.cc:642
 hparticletracksorter.cc:643
 hparticletracksorter.cc:644
 hparticletracksorter.cc:645
 hparticletracksorter.cc:646
 hparticletracksorter.cc:647
 hparticletracksorter.cc:648
 hparticletracksorter.cc:649
 hparticletracksorter.cc:650
 hparticletracksorter.cc:651
 hparticletracksorter.cc:652
 hparticletracksorter.cc:653
 hparticletracksorter.cc:654
 hparticletracksorter.cc:655
 hparticletracksorter.cc:656
 hparticletracksorter.cc:657
 hparticletracksorter.cc:658
 hparticletracksorter.cc:659
 hparticletracksorter.cc:660
 hparticletracksorter.cc:661
 hparticletracksorter.cc:662
 hparticletracksorter.cc:663
 hparticletracksorter.cc:664
 hparticletracksorter.cc:665
 hparticletracksorter.cc:666
 hparticletracksorter.cc:667
 hparticletracksorter.cc:668
 hparticletracksorter.cc:669
 hparticletracksorter.cc:670
 hparticletracksorter.cc:671
 hparticletracksorter.cc:672
 hparticletracksorter.cc:673
 hparticletracksorter.cc:674
 hparticletracksorter.cc:675
 hparticletracksorter.cc:676
 hparticletracksorter.cc:677
 hparticletracksorter.cc:678
 hparticletracksorter.cc:679
 hparticletracksorter.cc:680
 hparticletracksorter.cc:681
 hparticletracksorter.cc:682
 hparticletracksorter.cc:683
 hparticletracksorter.cc:684
 hparticletracksorter.cc:685
 hparticletracksorter.cc:686
 hparticletracksorter.cc:687
 hparticletracksorter.cc:688
 hparticletracksorter.cc:689
 hparticletracksorter.cc:690
 hparticletracksorter.cc:691
 hparticletracksorter.cc:692
 hparticletracksorter.cc:693
 hparticletracksorter.cc:694
 hparticletracksorter.cc:695
 hparticletracksorter.cc:696
 hparticletracksorter.cc:697
 hparticletracksorter.cc:698
 hparticletracksorter.cc:699
 hparticletracksorter.cc:700
 hparticletracksorter.cc:701
 hparticletracksorter.cc:702
 hparticletracksorter.cc:703
 hparticletracksorter.cc:704
 hparticletracksorter.cc:705
 hparticletracksorter.cc:706
 hparticletracksorter.cc:707
 hparticletracksorter.cc:708
 hparticletracksorter.cc:709
 hparticletracksorter.cc:710
 hparticletracksorter.cc:711
 hparticletracksorter.cc:712
 hparticletracksorter.cc:713
 hparticletracksorter.cc:714
 hparticletracksorter.cc:715
 hparticletracksorter.cc:716
 hparticletracksorter.cc:717
 hparticletracksorter.cc:718
 hparticletracksorter.cc:719
 hparticletracksorter.cc:720
 hparticletracksorter.cc:721
 hparticletracksorter.cc:722
 hparticletracksorter.cc:723
 hparticletracksorter.cc:724
 hparticletracksorter.cc:725
 hparticletracksorter.cc:726
 hparticletracksorter.cc:727
 hparticletracksorter.cc:728
 hparticletracksorter.cc:729
 hparticletracksorter.cc:730
 hparticletracksorter.cc:731
 hparticletracksorter.cc:732
 hparticletracksorter.cc:733
 hparticletracksorter.cc:734
 hparticletracksorter.cc:735
 hparticletracksorter.cc:736
 hparticletracksorter.cc:737
 hparticletracksorter.cc:738
 hparticletracksorter.cc:739
 hparticletracksorter.cc:740
 hparticletracksorter.cc:741
 hparticletracksorter.cc:742
 hparticletracksorter.cc:743
 hparticletracksorter.cc:744
 hparticletracksorter.cc:745
 hparticletracksorter.cc:746
 hparticletracksorter.cc:747
 hparticletracksorter.cc:748
 hparticletracksorter.cc:749
 hparticletracksorter.cc:750
 hparticletracksorter.cc:751
 hparticletracksorter.cc:752
 hparticletracksorter.cc:753
 hparticletracksorter.cc:754
 hparticletracksorter.cc:755
 hparticletracksorter.cc:756
 hparticletracksorter.cc:757
 hparticletracksorter.cc:758
 hparticletracksorter.cc:759
 hparticletracksorter.cc:760
 hparticletracksorter.cc:761
 hparticletracksorter.cc:762
 hparticletracksorter.cc:763
 hparticletracksorter.cc:764
 hparticletracksorter.cc:765
 hparticletracksorter.cc:766
 hparticletracksorter.cc:767
 hparticletracksorter.cc:768
 hparticletracksorter.cc:769
 hparticletracksorter.cc:770
 hparticletracksorter.cc:771
 hparticletracksorter.cc:772
 hparticletracksorter.cc:773
 hparticletracksorter.cc:774
 hparticletracksorter.cc:775
 hparticletracksorter.cc:776
 hparticletracksorter.cc:777
 hparticletracksorter.cc:778
 hparticletracksorter.cc:779
 hparticletracksorter.cc:780
 hparticletracksorter.cc:781
 hparticletracksorter.cc:782
 hparticletracksorter.cc:783
 hparticletracksorter.cc:784
 hparticletracksorter.cc:785
 hparticletracksorter.cc:786
 hparticletracksorter.cc:787
 hparticletracksorter.cc:788
 hparticletracksorter.cc:789
 hparticletracksorter.cc:790
 hparticletracksorter.cc:791
 hparticletracksorter.cc:792
 hparticletracksorter.cc:793
 hparticletracksorter.cc:794
 hparticletracksorter.cc:795
 hparticletracksorter.cc:796
 hparticletracksorter.cc:797
 hparticletracksorter.cc:798
 hparticletracksorter.cc:799
 hparticletracksorter.cc:800
 hparticletracksorter.cc:801
 hparticletracksorter.cc:802
 hparticletracksorter.cc:803
 hparticletracksorter.cc:804
 hparticletracksorter.cc:805
 hparticletracksorter.cc:806
 hparticletracksorter.cc:807
 hparticletracksorter.cc:808
 hparticletracksorter.cc:809
 hparticletracksorter.cc:810
 hparticletracksorter.cc:811
 hparticletracksorter.cc:812
 hparticletracksorter.cc:813
 hparticletracksorter.cc:814
 hparticletracksorter.cc:815
 hparticletracksorter.cc:816
 hparticletracksorter.cc:817
 hparticletracksorter.cc:818
 hparticletracksorter.cc:819
 hparticletracksorter.cc:820
 hparticletracksorter.cc:821
 hparticletracksorter.cc:822
 hparticletracksorter.cc:823
 hparticletracksorter.cc:824
 hparticletracksorter.cc:825
 hparticletracksorter.cc:826
 hparticletracksorter.cc:827
 hparticletracksorter.cc:828
 hparticletracksorter.cc:829
 hparticletracksorter.cc:830
 hparticletracksorter.cc:831
 hparticletracksorter.cc:832
 hparticletracksorter.cc:833
 hparticletracksorter.cc:834
 hparticletracksorter.cc:835
 hparticletracksorter.cc:836
 hparticletracksorter.cc:837
 hparticletracksorter.cc:838
 hparticletracksorter.cc:839
 hparticletracksorter.cc:840
 hparticletracksorter.cc:841
 hparticletracksorter.cc:842
 hparticletracksorter.cc:843
 hparticletracksorter.cc:844
 hparticletracksorter.cc:845
 hparticletracksorter.cc:846
 hparticletracksorter.cc:847
 hparticletracksorter.cc:848
 hparticletracksorter.cc:849
 hparticletracksorter.cc:850
 hparticletracksorter.cc:851
 hparticletracksorter.cc:852
 hparticletracksorter.cc:853
 hparticletracksorter.cc:854
 hparticletracksorter.cc:855
 hparticletracksorter.cc:856
 hparticletracksorter.cc:857
 hparticletracksorter.cc:858
 hparticletracksorter.cc:859
 hparticletracksorter.cc:860
 hparticletracksorter.cc:861
 hparticletracksorter.cc:862
 hparticletracksorter.cc:863
 hparticletracksorter.cc:864
 hparticletracksorter.cc:865
 hparticletracksorter.cc:866
 hparticletracksorter.cc:867
 hparticletracksorter.cc:868
 hparticletracksorter.cc:869
 hparticletracksorter.cc:870
 hparticletracksorter.cc:871
 hparticletracksorter.cc:872
 hparticletracksorter.cc:873
 hparticletracksorter.cc:874
 hparticletracksorter.cc:875
 hparticletracksorter.cc:876
 hparticletracksorter.cc:877
 hparticletracksorter.cc:878
 hparticletracksorter.cc:879
 hparticletracksorter.cc:880
 hparticletracksorter.cc:881
 hparticletracksorter.cc:882
 hparticletracksorter.cc:883
 hparticletracksorter.cc:884
 hparticletracksorter.cc:885
 hparticletracksorter.cc:886
 hparticletracksorter.cc:887
 hparticletracksorter.cc:888
 hparticletracksorter.cc:889
 hparticletracksorter.cc:890
 hparticletracksorter.cc:891
 hparticletracksorter.cc:892
 hparticletracksorter.cc:893
 hparticletracksorter.cc:894
 hparticletracksorter.cc:895
 hparticletracksorter.cc:896
 hparticletracksorter.cc:897
 hparticletracksorter.cc:898
 hparticletracksorter.cc:899
 hparticletracksorter.cc:900
 hparticletracksorter.cc:901
 hparticletracksorter.cc:902
 hparticletracksorter.cc:903
 hparticletracksorter.cc:904
 hparticletracksorter.cc:905
 hparticletracksorter.cc:906
 hparticletracksorter.cc:907
 hparticletracksorter.cc:908
 hparticletracksorter.cc:909
 hparticletracksorter.cc:910
 hparticletracksorter.cc:911
 hparticletracksorter.cc:912
 hparticletracksorter.cc:913
 hparticletracksorter.cc:914
 hparticletracksorter.cc:915
 hparticletracksorter.cc:916
 hparticletracksorter.cc:917
 hparticletracksorter.cc:918
 hparticletracksorter.cc:919
 hparticletracksorter.cc:920
 hparticletracksorter.cc:921
 hparticletracksorter.cc:922
 hparticletracksorter.cc:923
 hparticletracksorter.cc:924
 hparticletracksorter.cc:925
 hparticletracksorter.cc:926
 hparticletracksorter.cc:927
 hparticletracksorter.cc:928
 hparticletracksorter.cc:929
 hparticletracksorter.cc:930
 hparticletracksorter.cc:931
 hparticletracksorter.cc:932
 hparticletracksorter.cc:933
 hparticletracksorter.cc:934
 hparticletracksorter.cc:935
 hparticletracksorter.cc:936
 hparticletracksorter.cc:937
 hparticletracksorter.cc:938
 hparticletracksorter.cc:939
 hparticletracksorter.cc:940
 hparticletracksorter.cc:941
 hparticletracksorter.cc:942
 hparticletracksorter.cc:943
 hparticletracksorter.cc:944
 hparticletracksorter.cc:945
 hparticletracksorter.cc:946
 hparticletracksorter.cc:947
 hparticletracksorter.cc:948
 hparticletracksorter.cc:949
 hparticletracksorter.cc:950
 hparticletracksorter.cc:951
 hparticletracksorter.cc:952
 hparticletracksorter.cc:953
 hparticletracksorter.cc:954
 hparticletracksorter.cc:955
 hparticletracksorter.cc:956
 hparticletracksorter.cc:957
 hparticletracksorter.cc:958
 hparticletracksorter.cc:959
 hparticletracksorter.cc:960
 hparticletracksorter.cc:961
 hparticletracksorter.cc:962
 hparticletracksorter.cc:963
 hparticletracksorter.cc:964
 hparticletracksorter.cc:965
 hparticletracksorter.cc:966
 hparticletracksorter.cc:967
 hparticletracksorter.cc:968
 hparticletracksorter.cc:969
 hparticletracksorter.cc:970
 hparticletracksorter.cc:971
 hparticletracksorter.cc:972
 hparticletracksorter.cc:973
 hparticletracksorter.cc:974
 hparticletracksorter.cc:975
 hparticletracksorter.cc:976
 hparticletracksorter.cc:977
 hparticletracksorter.cc:978
 hparticletracksorter.cc:979
 hparticletracksorter.cc:980
 hparticletracksorter.cc:981
 hparticletracksorter.cc:982
 hparticletracksorter.cc:983
 hparticletracksorter.cc:984
 hparticletracksorter.cc:985
 hparticletracksorter.cc:986
 hparticletracksorter.cc:987
 hparticletracksorter.cc:988
 hparticletracksorter.cc:989
 hparticletracksorter.cc:990
 hparticletracksorter.cc:991
 hparticletracksorter.cc:992
 hparticletracksorter.cc:993
 hparticletracksorter.cc:994
 hparticletracksorter.cc:995
 hparticletracksorter.cc:996
 hparticletracksorter.cc:997
 hparticletracksorter.cc:998
 hparticletracksorter.cc:999
 hparticletracksorter.cc:1000
 hparticletracksorter.cc:1001
 hparticletracksorter.cc:1002
 hparticletracksorter.cc:1003
 hparticletracksorter.cc:1004
 hparticletracksorter.cc:1005
 hparticletracksorter.cc:1006
 hparticletracksorter.cc:1007
 hparticletracksorter.cc:1008
 hparticletracksorter.cc:1009
 hparticletracksorter.cc:1010
 hparticletracksorter.cc:1011
 hparticletracksorter.cc:1012
 hparticletracksorter.cc:1013
 hparticletracksorter.cc:1014
 hparticletracksorter.cc:1015
 hparticletracksorter.cc:1016
 hparticletracksorter.cc:1017
 hparticletracksorter.cc:1018
 hparticletracksorter.cc:1019
 hparticletracksorter.cc:1020
 hparticletracksorter.cc:1021
 hparticletracksorter.cc:1022
 hparticletracksorter.cc:1023
 hparticletracksorter.cc:1024
 hparticletracksorter.cc:1025
 hparticletracksorter.cc:1026
 hparticletracksorter.cc:1027
 hparticletracksorter.cc:1028
 hparticletracksorter.cc:1029
 hparticletracksorter.cc:1030
 hparticletracksorter.cc:1031
 hparticletracksorter.cc:1032
 hparticletracksorter.cc:1033
 hparticletracksorter.cc:1034
 hparticletracksorter.cc:1035
 hparticletracksorter.cc:1036
 hparticletracksorter.cc:1037
 hparticletracksorter.cc:1038
 hparticletracksorter.cc:1039
 hparticletracksorter.cc:1040
 hparticletracksorter.cc:1041
 hparticletracksorter.cc:1042
 hparticletracksorter.cc:1043
 hparticletracksorter.cc:1044
 hparticletracksorter.cc:1045
 hparticletracksorter.cc:1046
 hparticletracksorter.cc:1047
 hparticletracksorter.cc:1048
 hparticletracksorter.cc:1049
 hparticletracksorter.cc:1050
 hparticletracksorter.cc:1051
 hparticletracksorter.cc:1052
 hparticletracksorter.cc:1053
 hparticletracksorter.cc:1054
 hparticletracksorter.cc:1055
 hparticletracksorter.cc:1056
 hparticletracksorter.cc:1057
 hparticletracksorter.cc:1058
 hparticletracksorter.cc:1059
 hparticletracksorter.cc:1060
 hparticletracksorter.cc:1061
 hparticletracksorter.cc:1062
 hparticletracksorter.cc:1063
 hparticletracksorter.cc:1064
 hparticletracksorter.cc:1065
 hparticletracksorter.cc:1066
 hparticletracksorter.cc:1067
 hparticletracksorter.cc:1068
 hparticletracksorter.cc:1069
 hparticletracksorter.cc:1070
 hparticletracksorter.cc:1071
 hparticletracksorter.cc:1072
 hparticletracksorter.cc:1073
 hparticletracksorter.cc:1074
 hparticletracksorter.cc:1075
 hparticletracksorter.cc:1076
 hparticletracksorter.cc:1077
 hparticletracksorter.cc:1078
 hparticletracksorter.cc:1079
 hparticletracksorter.cc:1080
 hparticletracksorter.cc:1081
 hparticletracksorter.cc:1082
 hparticletracksorter.cc:1083
 hparticletracksorter.cc:1084
 hparticletracksorter.cc:1085
 hparticletracksorter.cc:1086
 hparticletracksorter.cc:1087
 hparticletracksorter.cc:1088
 hparticletracksorter.cc:1089
 hparticletracksorter.cc:1090
 hparticletracksorter.cc:1091
 hparticletracksorter.cc:1092
 hparticletracksorter.cc:1093
 hparticletracksorter.cc:1094
 hparticletracksorter.cc:1095
 hparticletracksorter.cc:1096
 hparticletracksorter.cc:1097
 hparticletracksorter.cc:1098
 hparticletracksorter.cc:1099
 hparticletracksorter.cc:1100
 hparticletracksorter.cc:1101
 hparticletracksorter.cc:1102
 hparticletracksorter.cc:1103
 hparticletracksorter.cc:1104
 hparticletracksorter.cc:1105
 hparticletracksorter.cc:1106
 hparticletracksorter.cc:1107
 hparticletracksorter.cc:1108
 hparticletracksorter.cc:1109
 hparticletracksorter.cc:1110
 hparticletracksorter.cc:1111
 hparticletracksorter.cc:1112
 hparticletracksorter.cc:1113
 hparticletracksorter.cc:1114
 hparticletracksorter.cc:1115
 hparticletracksorter.cc:1116
 hparticletracksorter.cc:1117
 hparticletracksorter.cc:1118
 hparticletracksorter.cc:1119
 hparticletracksorter.cc:1120
 hparticletracksorter.cc:1121
 hparticletracksorter.cc:1122
 hparticletracksorter.cc:1123
 hparticletracksorter.cc:1124
 hparticletracksorter.cc:1125
 hparticletracksorter.cc:1126
 hparticletracksorter.cc:1127
 hparticletracksorter.cc:1128
 hparticletracksorter.cc:1129
 hparticletracksorter.cc:1130
 hparticletracksorter.cc:1131
 hparticletracksorter.cc:1132
 hparticletracksorter.cc:1133
 hparticletracksorter.cc:1134
 hparticletracksorter.cc:1135
 hparticletracksorter.cc:1136
 hparticletracksorter.cc:1137
 hparticletracksorter.cc:1138
 hparticletracksorter.cc:1139
 hparticletracksorter.cc:1140
 hparticletracksorter.cc:1141
 hparticletracksorter.cc:1142
 hparticletracksorter.cc:1143
 hparticletracksorter.cc:1144
 hparticletracksorter.cc:1145
 hparticletracksorter.cc:1146
 hparticletracksorter.cc:1147
 hparticletracksorter.cc:1148
 hparticletracksorter.cc:1149
 hparticletracksorter.cc:1150
 hparticletracksorter.cc:1151
 hparticletracksorter.cc:1152
 hparticletracksorter.cc:1153
 hparticletracksorter.cc:1154
 hparticletracksorter.cc:1155
 hparticletracksorter.cc:1156
 hparticletracksorter.cc:1157
 hparticletracksorter.cc:1158
 hparticletracksorter.cc:1159
 hparticletracksorter.cc:1160
 hparticletracksorter.cc:1161
 hparticletracksorter.cc:1162
 hparticletracksorter.cc:1163
 hparticletracksorter.cc:1164
 hparticletracksorter.cc:1165
 hparticletracksorter.cc:1166
 hparticletracksorter.cc:1167
 hparticletracksorter.cc:1168
 hparticletracksorter.cc:1169
 hparticletracksorter.cc:1170
 hparticletracksorter.cc:1171
 hparticletracksorter.cc:1172
 hparticletracksorter.cc:1173
 hparticletracksorter.cc:1174
 hparticletracksorter.cc:1175
 hparticletracksorter.cc:1176
 hparticletracksorter.cc:1177
 hparticletracksorter.cc:1178
 hparticletracksorter.cc:1179
 hparticletracksorter.cc:1180
 hparticletracksorter.cc:1181
 hparticletracksorter.cc:1182
 hparticletracksorter.cc:1183
 hparticletracksorter.cc:1184
 hparticletracksorter.cc:1185
 hparticletracksorter.cc:1186
 hparticletracksorter.cc:1187
 hparticletracksorter.cc:1188
 hparticletracksorter.cc:1189
 hparticletracksorter.cc:1190
 hparticletracksorter.cc:1191
 hparticletracksorter.cc:1192
 hparticletracksorter.cc:1193
 hparticletracksorter.cc:1194
 hparticletracksorter.cc:1195
 hparticletracksorter.cc:1196
 hparticletracksorter.cc:1197
 hparticletracksorter.cc:1198
 hparticletracksorter.cc:1199
 hparticletracksorter.cc:1200
 hparticletracksorter.cc:1201
 hparticletracksorter.cc:1202
 hparticletracksorter.cc:1203
 hparticletracksorter.cc:1204
 hparticletracksorter.cc:1205
 hparticletracksorter.cc:1206
 hparticletracksorter.cc:1207
 hparticletracksorter.cc:1208
 hparticletracksorter.cc:1209
 hparticletracksorter.cc:1210
 hparticletracksorter.cc:1211
 hparticletracksorter.cc:1212
 hparticletracksorter.cc:1213
 hparticletracksorter.cc:1214
 hparticletracksorter.cc:1215
 hparticletracksorter.cc:1216
 hparticletracksorter.cc:1217
 hparticletracksorter.cc:1218
 hparticletracksorter.cc:1219
 hparticletracksorter.cc:1220
 hparticletracksorter.cc:1221
 hparticletracksorter.cc:1222
 hparticletracksorter.cc:1223
 hparticletracksorter.cc:1224
 hparticletracksorter.cc:1225
 hparticletracksorter.cc:1226
 hparticletracksorter.cc:1227
 hparticletracksorter.cc:1228
 hparticletracksorter.cc:1229
 hparticletracksorter.cc:1230
 hparticletracksorter.cc:1231
 hparticletracksorter.cc:1232
 hparticletracksorter.cc:1233
 hparticletracksorter.cc:1234
 hparticletracksorter.cc:1235
 hparticletracksorter.cc:1236
 hparticletracksorter.cc:1237
 hparticletracksorter.cc:1238
 hparticletracksorter.cc:1239
 hparticletracksorter.cc:1240
 hparticletracksorter.cc:1241
 hparticletracksorter.cc:1242
 hparticletracksorter.cc:1243
 hparticletracksorter.cc:1244
 hparticletracksorter.cc:1245
 hparticletracksorter.cc:1246
 hparticletracksorter.cc:1247
 hparticletracksorter.cc:1248
 hparticletracksorter.cc:1249
 hparticletracksorter.cc:1250
 hparticletracksorter.cc:1251
 hparticletracksorter.cc:1252
 hparticletracksorter.cc:1253
 hparticletracksorter.cc:1254
 hparticletracksorter.cc:1255
 hparticletracksorter.cc:1256
 hparticletracksorter.cc:1257
 hparticletracksorter.cc:1258
 hparticletracksorter.cc:1259
 hparticletracksorter.cc:1260
 hparticletracksorter.cc:1261
 hparticletracksorter.cc:1262
 hparticletracksorter.cc:1263
 hparticletracksorter.cc:1264
 hparticletracksorter.cc:1265
 hparticletracksorter.cc:1266
 hparticletracksorter.cc:1267
 hparticletracksorter.cc:1268
 hparticletracksorter.cc:1269
 hparticletracksorter.cc:1270
 hparticletracksorter.cc:1271
 hparticletracksorter.cc:1272
 hparticletracksorter.cc:1273
 hparticletracksorter.cc:1274
 hparticletracksorter.cc:1275
 hparticletracksorter.cc:1276
 hparticletracksorter.cc:1277
 hparticletracksorter.cc:1278
 hparticletracksorter.cc:1279
 hparticletracksorter.cc:1280
 hparticletracksorter.cc:1281
 hparticletracksorter.cc:1282
 hparticletracksorter.cc:1283
 hparticletracksorter.cc:1284
 hparticletracksorter.cc:1285
 hparticletracksorter.cc:1286
 hparticletracksorter.cc:1287
 hparticletracksorter.cc:1288
 hparticletracksorter.cc:1289
 hparticletracksorter.cc:1290
 hparticletracksorter.cc:1291
 hparticletracksorter.cc:1292
 hparticletracksorter.cc:1293
 hparticletracksorter.cc:1294
 hparticletracksorter.cc:1295
 hparticletracksorter.cc:1296
 hparticletracksorter.cc:1297
 hparticletracksorter.cc:1298
 hparticletracksorter.cc:1299
 hparticletracksorter.cc:1300
 hparticletracksorter.cc:1301
 hparticletracksorter.cc:1302
 hparticletracksorter.cc:1303
 hparticletracksorter.cc:1304
 hparticletracksorter.cc:1305
 hparticletracksorter.cc:1306
 hparticletracksorter.cc:1307
 hparticletracksorter.cc:1308
 hparticletracksorter.cc:1309
 hparticletracksorter.cc:1310
 hparticletracksorter.cc:1311
 hparticletracksorter.cc:1312
 hparticletracksorter.cc:1313
 hparticletracksorter.cc:1314
 hparticletracksorter.cc:1315
 hparticletracksorter.cc:1316
 hparticletracksorter.cc:1317
 hparticletracksorter.cc:1318
 hparticletracksorter.cc:1319
 hparticletracksorter.cc:1320
 hparticletracksorter.cc:1321
 hparticletracksorter.cc:1322
 hparticletracksorter.cc:1323
 hparticletracksorter.cc:1324
 hparticletracksorter.cc:1325
 hparticletracksorter.cc:1326
 hparticletracksorter.cc:1327
 hparticletracksorter.cc:1328
 hparticletracksorter.cc:1329
 hparticletracksorter.cc:1330
 hparticletracksorter.cc:1331
 hparticletracksorter.cc:1332
 hparticletracksorter.cc:1333
 hparticletracksorter.cc:1334
 hparticletracksorter.cc:1335
 hparticletracksorter.cc:1336
 hparticletracksorter.cc:1337
 hparticletracksorter.cc:1338
 hparticletracksorter.cc:1339
 hparticletracksorter.cc:1340
 hparticletracksorter.cc:1341
 hparticletracksorter.cc:1342
 hparticletracksorter.cc:1343
 hparticletracksorter.cc:1344
 hparticletracksorter.cc:1345
 hparticletracksorter.cc:1346
 hparticletracksorter.cc:1347
 hparticletracksorter.cc:1348
 hparticletracksorter.cc:1349
 hparticletracksorter.cc:1350
 hparticletracksorter.cc:1351
 hparticletracksorter.cc:1352
 hparticletracksorter.cc:1353
 hparticletracksorter.cc:1354
 hparticletracksorter.cc:1355
 hparticletracksorter.cc:1356
 hparticletracksorter.cc:1357
 hparticletracksorter.cc:1358
 hparticletracksorter.cc:1359
 hparticletracksorter.cc:1360
 hparticletracksorter.cc:1361
 hparticletracksorter.cc:1362
 hparticletracksorter.cc:1363
 hparticletracksorter.cc:1364
 hparticletracksorter.cc:1365
 hparticletracksorter.cc:1366
 hparticletracksorter.cc:1367
 hparticletracksorter.cc:1368
 hparticletracksorter.cc:1369
 hparticletracksorter.cc:1370
 hparticletracksorter.cc:1371
 hparticletracksorter.cc:1372
 hparticletracksorter.cc:1373
 hparticletracksorter.cc:1374
 hparticletracksorter.cc:1375
 hparticletracksorter.cc:1376
 hparticletracksorter.cc:1377
 hparticletracksorter.cc:1378
 hparticletracksorter.cc:1379
 hparticletracksorter.cc:1380
 hparticletracksorter.cc:1381
 hparticletracksorter.cc:1382
 hparticletracksorter.cc:1383
 hparticletracksorter.cc:1384
 hparticletracksorter.cc:1385
 hparticletracksorter.cc:1386
 hparticletracksorter.cc:1387
 hparticletracksorter.cc:1388
 hparticletracksorter.cc:1389
 hparticletracksorter.cc:1390
 hparticletracksorter.cc:1391
 hparticletracksorter.cc:1392
 hparticletracksorter.cc:1393
 hparticletracksorter.cc:1394
 hparticletracksorter.cc:1395
 hparticletracksorter.cc:1396
 hparticletracksorter.cc:1397
 hparticletracksorter.cc:1398
 hparticletracksorter.cc:1399
 hparticletracksorter.cc:1400
 hparticletracksorter.cc:1401
 hparticletracksorter.cc:1402
 hparticletracksorter.cc:1403
 hparticletracksorter.cc:1404
 hparticletracksorter.cc:1405
 hparticletracksorter.cc:1406
 hparticletracksorter.cc:1407
 hparticletracksorter.cc:1408
 hparticletracksorter.cc:1409
 hparticletracksorter.cc:1410
 hparticletracksorter.cc:1411
 hparticletracksorter.cc:1412
 hparticletracksorter.cc:1413
 hparticletracksorter.cc:1414
 hparticletracksorter.cc:1415
 hparticletracksorter.cc:1416
 hparticletracksorter.cc:1417
 hparticletracksorter.cc:1418
 hparticletracksorter.cc:1419
 hparticletracksorter.cc:1420
 hparticletracksorter.cc:1421
 hparticletracksorter.cc:1422
 hparticletracksorter.cc:1423
 hparticletracksorter.cc:1424
 hparticletracksorter.cc:1425
 hparticletracksorter.cc:1426
 hparticletracksorter.cc:1427
 hparticletracksorter.cc:1428
 hparticletracksorter.cc:1429
 hparticletracksorter.cc:1430
 hparticletracksorter.cc:1431
 hparticletracksorter.cc:1432
 hparticletracksorter.cc:1433
 hparticletracksorter.cc:1434
 hparticletracksorter.cc:1435
 hparticletracksorter.cc:1436
 hparticletracksorter.cc:1437
 hparticletracksorter.cc:1438
 hparticletracksorter.cc:1439
 hparticletracksorter.cc:1440
 hparticletracksorter.cc:1441
 hparticletracksorter.cc:1442
 hparticletracksorter.cc:1443
 hparticletracksorter.cc:1444
 hparticletracksorter.cc:1445
 hparticletracksorter.cc:1446
 hparticletracksorter.cc:1447
 hparticletracksorter.cc:1448
 hparticletracksorter.cc:1449
 hparticletracksorter.cc:1450
 hparticletracksorter.cc:1451
 hparticletracksorter.cc:1452
 hparticletracksorter.cc:1453
 hparticletracksorter.cc:1454
 hparticletracksorter.cc:1455
 hparticletracksorter.cc:1456
 hparticletracksorter.cc:1457
 hparticletracksorter.cc:1458
 hparticletracksorter.cc:1459
 hparticletracksorter.cc:1460
 hparticletracksorter.cc:1461
 hparticletracksorter.cc:1462
 hparticletracksorter.cc:1463
 hparticletracksorter.cc:1464
 hparticletracksorter.cc:1465
 hparticletracksorter.cc:1466
 hparticletracksorter.cc:1467
 hparticletracksorter.cc:1468
 hparticletracksorter.cc:1469
 hparticletracksorter.cc:1470
 hparticletracksorter.cc:1471
 hparticletracksorter.cc:1472
 hparticletracksorter.cc:1473
 hparticletracksorter.cc:1474
 hparticletracksorter.cc:1475
 hparticletracksorter.cc:1476
 hparticletracksorter.cc:1477
 hparticletracksorter.cc:1478
 hparticletracksorter.cc:1479
 hparticletracksorter.cc:1480
 hparticletracksorter.cc:1481
 hparticletracksorter.cc:1482
 hparticletracksorter.cc:1483
 hparticletracksorter.cc:1484
 hparticletracksorter.cc:1485
 hparticletracksorter.cc:1486
 hparticletracksorter.cc:1487
 hparticletracksorter.cc:1488
 hparticletracksorter.cc:1489
 hparticletracksorter.cc:1490
 hparticletracksorter.cc:1491
 hparticletracksorter.cc:1492
 hparticletracksorter.cc:1493
 hparticletracksorter.cc:1494
 hparticletracksorter.cc:1495
 hparticletracksorter.cc:1496
 hparticletracksorter.cc:1497
 hparticletracksorter.cc:1498
 hparticletracksorter.cc:1499
 hparticletracksorter.cc:1500
 hparticletracksorter.cc:1501
 hparticletracksorter.cc:1502
 hparticletracksorter.cc:1503
 hparticletracksorter.cc:1504
 hparticletracksorter.cc:1505
 hparticletracksorter.cc:1506
 hparticletracksorter.cc:1507
 hparticletracksorter.cc:1508
 hparticletracksorter.cc:1509
 hparticletracksorter.cc:1510
 hparticletracksorter.cc:1511
 hparticletracksorter.cc:1512
 hparticletracksorter.cc:1513
 hparticletracksorter.cc:1514
 hparticletracksorter.cc:1515
 hparticletracksorter.cc:1516
 hparticletracksorter.cc:1517
 hparticletracksorter.cc:1518
 hparticletracksorter.cc:1519
 hparticletracksorter.cc:1520
 hparticletracksorter.cc:1521
 hparticletracksorter.cc:1522
 hparticletracksorter.cc:1523
 hparticletracksorter.cc:1524
 hparticletracksorter.cc:1525
 hparticletracksorter.cc:1526
 hparticletracksorter.cc:1527
 hparticletracksorter.cc:1528
 hparticletracksorter.cc:1529
 hparticletracksorter.cc:1530
 hparticletracksorter.cc:1531
 hparticletracksorter.cc:1532
 hparticletracksorter.cc:1533
 hparticletracksorter.cc:1534
 hparticletracksorter.cc:1535
 hparticletracksorter.cc:1536
 hparticletracksorter.cc:1537
 hparticletracksorter.cc:1538
 hparticletracksorter.cc:1539
 hparticletracksorter.cc:1540
 hparticletracksorter.cc:1541
 hparticletracksorter.cc:1542
 hparticletracksorter.cc:1543
 hparticletracksorter.cc:1544
 hparticletracksorter.cc:1545
 hparticletracksorter.cc:1546
 hparticletracksorter.cc:1547
 hparticletracksorter.cc:1548
 hparticletracksorter.cc:1549
 hparticletracksorter.cc:1550
 hparticletracksorter.cc:1551
 hparticletracksorter.cc:1552
 hparticletracksorter.cc:1553
 hparticletracksorter.cc:1554
 hparticletracksorter.cc:1555
 hparticletracksorter.cc:1556
 hparticletracksorter.cc:1557
 hparticletracksorter.cc:1558
 hparticletracksorter.cc:1559
 hparticletracksorter.cc:1560
 hparticletracksorter.cc:1561
 hparticletracksorter.cc:1562
 hparticletracksorter.cc:1563
 hparticletracksorter.cc:1564
 hparticletracksorter.cc:1565
 hparticletracksorter.cc:1566
 hparticletracksorter.cc:1567
 hparticletracksorter.cc:1568
 hparticletracksorter.cc:1569
 hparticletracksorter.cc:1570
 hparticletracksorter.cc:1571
 hparticletracksorter.cc:1572
 hparticletracksorter.cc:1573
 hparticletracksorter.cc:1574
 hparticletracksorter.cc:1575
 hparticletracksorter.cc:1576
 hparticletracksorter.cc:1577
 hparticletracksorter.cc:1578
 hparticletracksorter.cc:1579
 hparticletracksorter.cc:1580
 hparticletracksorter.cc:1581
 hparticletracksorter.cc:1582
 hparticletracksorter.cc:1583
 hparticletracksorter.cc:1584
 hparticletracksorter.cc:1585
 hparticletracksorter.cc:1586
 hparticletracksorter.cc:1587
 hparticletracksorter.cc:1588
 hparticletracksorter.cc:1589
 hparticletracksorter.cc:1590
 hparticletracksorter.cc:1591
 hparticletracksorter.cc:1592
 hparticletracksorter.cc:1593
 hparticletracksorter.cc:1594
 hparticletracksorter.cc:1595
 hparticletracksorter.cc:1596
 hparticletracksorter.cc:1597
 hparticletracksorter.cc:1598
 hparticletracksorter.cc:1599
 hparticletracksorter.cc:1600
 hparticletracksorter.cc:1601
 hparticletracksorter.cc:1602
 hparticletracksorter.cc:1603
 hparticletracksorter.cc:1604
 hparticletracksorter.cc:1605
 hparticletracksorter.cc:1606
 hparticletracksorter.cc:1607
 hparticletracksorter.cc:1608
 hparticletracksorter.cc:1609
 hparticletracksorter.cc:1610
 hparticletracksorter.cc:1611
 hparticletracksorter.cc:1612
 hparticletracksorter.cc:1613
 hparticletracksorter.cc:1614
 hparticletracksorter.cc:1615
 hparticletracksorter.cc:1616
 hparticletracksorter.cc:1617
 hparticletracksorter.cc:1618
 hparticletracksorter.cc:1619
 hparticletracksorter.cc:1620
 hparticletracksorter.cc:1621
 hparticletracksorter.cc:1622
 hparticletracksorter.cc:1623
 hparticletracksorter.cc:1624
 hparticletracksorter.cc:1625
 hparticletracksorter.cc:1626
 hparticletracksorter.cc:1627
 hparticletracksorter.cc:1628
 hparticletracksorter.cc:1629
 hparticletracksorter.cc:1630
 hparticletracksorter.cc:1631
 hparticletracksorter.cc:1632
 hparticletracksorter.cc:1633
 hparticletracksorter.cc:1634
 hparticletracksorter.cc:1635
 hparticletracksorter.cc:1636
 hparticletracksorter.cc:1637
 hparticletracksorter.cc:1638
 hparticletracksorter.cc:1639
 hparticletracksorter.cc:1640
 hparticletracksorter.cc:1641
 hparticletracksorter.cc:1642
 hparticletracksorter.cc:1643
 hparticletracksorter.cc:1644
 hparticletracksorter.cc:1645
 hparticletracksorter.cc:1646
 hparticletracksorter.cc:1647
 hparticletracksorter.cc:1648
 hparticletracksorter.cc:1649
 hparticletracksorter.cc:1650
 hparticletracksorter.cc:1651
 hparticletracksorter.cc:1652
 hparticletracksorter.cc:1653
 hparticletracksorter.cc:1654
 hparticletracksorter.cc:1655
 hparticletracksorter.cc:1656
 hparticletracksorter.cc:1657
 hparticletracksorter.cc:1658
 hparticletracksorter.cc:1659
 hparticletracksorter.cc:1660
 hparticletracksorter.cc:1661
 hparticletracksorter.cc:1662
 hparticletracksorter.cc:1663
 hparticletracksorter.cc:1664
 hparticletracksorter.cc:1665
 hparticletracksorter.cc:1666
 hparticletracksorter.cc:1667
 hparticletracksorter.cc:1668
 hparticletracksorter.cc:1669
 hparticletracksorter.cc:1670
 hparticletracksorter.cc:1671
 hparticletracksorter.cc:1672
 hparticletracksorter.cc:1673
 hparticletracksorter.cc:1674
 hparticletracksorter.cc:1675
 hparticletracksorter.cc:1676
 hparticletracksorter.cc:1677
 hparticletracksorter.cc:1678
 hparticletracksorter.cc:1679
 hparticletracksorter.cc:1680
 hparticletracksorter.cc:1681
 hparticletracksorter.cc:1682
 hparticletracksorter.cc:1683
 hparticletracksorter.cc:1684
 hparticletracksorter.cc:1685
 hparticletracksorter.cc:1686
 hparticletracksorter.cc:1687
 hparticletracksorter.cc:1688
 hparticletracksorter.cc:1689
 hparticletracksorter.cc:1690
 hparticletracksorter.cc:1691
 hparticletracksorter.cc:1692
 hparticletracksorter.cc:1693
 hparticletracksorter.cc:1694
 hparticletracksorter.cc:1695
 hparticletracksorter.cc:1696
 hparticletracksorter.cc:1697
 hparticletracksorter.cc:1698
 hparticletracksorter.cc:1699
 hparticletracksorter.cc:1700
 hparticletracksorter.cc:1701
 hparticletracksorter.cc:1702
 hparticletracksorter.cc:1703
 hparticletracksorter.cc:1704