ROOT logo
//*-- AUTHOR : Vladimir Pechenov
//*-- Modified : 02/02/07

//_HADES_CLASS_DESCRIPTION
////////////////////////////////////////////////////////////////
//
// HMdcGeantSeg
//
//   Class collect and keep GEANT and mdc digitized information about
//   one inner or outer mdc segment
//
//
// HMdcGeantTrack
//
//   Class collect and keep GEANT and mdc digitized information about
//   one track. Keep also array of HMdcGeantSeg objects for this track.
//
//
// HMdcGeantEvent
//
//   Class fill and keep array of HMdcGeantTrack objects for one event.
//
//
// HMdcTrackInfSim
//
//   Class fill HMdcClusSim, HMdcSegSim and HMdcHitSim objects
//   by information from HMdcGeantTrack and HMdcGeantSeg
//
//
// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
// Meaning of HMdcGeantTrack::trackStatus, HMdcGeantSeg::segmentStatus
// and HMdcClusSim::trackStatus[tr] bits.
//
// Bits in HMdcGeantTrack::trackStatus:
// ------------------------------------
//  Bit 8, GEANT_BUG          : =0 - GEANT bug was found in this track
//  Bit 7, NO_META_HIT        : =0 - track don't reach META
//  Bit 6, NOT_RECONSTRUCTABLE: =0 - not reconstructable track. It can be if
//                                   track not present in mdc,
//                                   track cross more then one sector,
//                                   numder of segments <1 or >2,
//                                   upstream direction of at list one segment,
//                                   track doesn't cross inner or outer segment,
//                                   num. of segments with digi.wires <1 or >2,
//                                   num.digi.wires in track < 5*num.segments,
//                                   at list one seg. has <3 layers or <5 wires
//  Bits 1-5 are eq.1
//
// Bits in HMdcGeantSeg::segmentStatus trackStatus bits plus:
// ----------------------------------------------------------
//  Bit 5, SEG_NOT_RECONS_ABLE: =0 - segment has <3 layers or <5 digi.wires
//  Bit 4, HIT_NOT_RECONS_ABLE: =0 - at list one mdc in segment has
//                                   <3 fired layers or <5 wires
//
// Bits in HMdcClusSim::trackStatus[tr] are segmentStatus bits plus:
// -----------------------------------------------------------------
//  Bit 3, CLFLEVEL_TOO_HIGH  : =0 - segment can't be found due to high level
//                                   of cluster finder
//  Bit 2, FAKE_CONTRIBUTION  : =0 - fake contribution in cluster
//  Bit 1, SEGS_MATH_OK       : =0 - no inner-outer segments matching
////////////////////////////////////////////////////////////////

using namespace std;
#include "hmdcgeanttrack.h"
#include "hades.h"
#include "hevent.h"
//#include "hcategory.h"
#include "hmatrixcategory.h"
#include "hiterator.h"
#include "hmdcsizescells.h"
#include "hgeantkine.h"
#include "hgeantmdc.h"
#include "hmdccal1sim.h"
#include "hmdcdef.h"
#include "hphysicsconstants.h"
#include "hspectrometer.h"
#include "hmdclistcells.h"
#include "hmdcclussim.h"
#include "hmdcsegsim.h"
#include "hmdchitsim.h"
#include "hmdcgetcontainers.h"
#include "hmdcwiredata.h"
#include <iostream>
#include <iomanip>
#include <stdlib.h>

ClassImp(HMdcGeantSeg)

ClassImp(HMdcGeantTrack)

ClassImp(HMdcGeantEvent)


#define DEFAULT_SET         255
#define SEGS_MATH_OK          1
#define GEANT_BUG           127
#define NO_META_HIT         191
#define NOT_RECONSTRUCTABLE 223
#define SEG_NOT_RECONS_ABLE 239
#define HIT_NOT_RECONS_ABLE 247
#define CLFLEVEL_TOO_HIGH   251
#define FAKE_CONTRIBUTION   253

HMdcGeantSeg::HMdcGeantSeg(Short_t ns) {
  clear(ns);
}

void HMdcGeantSeg::clear(Short_t ns) {
  nHitsMdc[0]    = 0;
  nHitsMdc[1]    = 0;
  segNumber      = ns;
  sec            = -1;
  mod            = -1;
  pClusBest      = NULL;
  pClusBestCh[0] = NULL;
  pClusBestCh[1] = NULL;
  pSegBest       = NULL;
  trackNumber    = -1;
  userFlag       = 0;
  for(Int_t m=0;m<2;m++) for(Int_t l=0;l<7;l++) geantLay[m][l]=0;
  HMdcList12GroupCells::clear();
}

void HMdcGeantSeg::addFirstHit(HGeantMdc* pGeantMdc,Bool_t* mdcSecSetup,Short_t ns) {
  for(Int_t m=0;m<2;m++) for(Int_t l=0;l<7;l++) geantLay[m][l]=0;
  HMdcList12GroupCells::clear();
  mod        = pGeantMdc->getModule();
  ioseg      = mod >> 1;
  Int_t modI = mod & 1;
  Int_t lay  = pGeantMdc->getLayer();
  sec        = pGeantMdc->getSector();
  segNumber  = ns;
  nHitsMdc[0]=nHitsMdc[1]=0;
  if(lay!=6) nHitsMdc[modI]++;
  geantLay[modI][lay] = pGeantMdc;
  direction = dirHit(pGeantMdc);
  if(mdcSecSetup[ioseg*2] && mdcSecSetup[ioseg*2+1]) mod = -2;
}

Char_t HMdcGeantSeg::dirHit(HGeantMdc* pGeantMdc) {
  Float_t ath, aph;
  pGeantMdc->getIncidence(ath,aph);
  return dirTheta(ath);
}

Bool_t HMdcGeantSeg::addHit(HGeantMdc* pGeantMdc) {
  if(pGeantMdc->getSector() != sec) return kFALSE;
  Char_t mod  = pGeantMdc->getModule();
  if(mod>>1 != ioseg) return kFALSE;
  if(direction != dirHit(pGeantMdc)) return kFALSE;
  Int_t modI = mod & 1;
  Int_t lay  = pGeantMdc->getLayer();
  if(geantLay[modI][lay]) return kFALSE;
  geantLay[modI][lay] = pGeantMdc;
  if(lay!=6) nHitsMdc[modI]++;
  return kTRUE;
}

Int_t HMdcGeantSeg::getModIOfGeantTrack(void) const {
  // return -2 if track present in both mdc of segment
  // otherwise return mdc index (0 or 1)
  if(nHitsMdc[0]>0 && nHitsMdc[1]>0) return -2;
  if(nHitsMdc[0]>0) return 0;
  if(nHitsMdc[1]>0) return 1;
  return -500;
}

void HMdcGeantSeg::setAnotherHit(HGeantMdc* pGeantMdc) {
  Int_t modI = pGeantMdc->getModule() & 1;  // =0 or 1 !
  Int_t lay  = pGeantMdc->getLayer();
  if(geantLay[modI][lay]) geantLay[modI][lay] = pGeantMdc;
}


Int_t HMdcGeantSeg::getFirstLayer12(void) const {
  if(nHitsMdc[0]) {
    for(Int_t l=0;l<6;l++) if(geantLay[0][l]) return l;
  } else for(Int_t l=0;l<6;l++) if(geantLay[1][l]) return l+6;
  return -1;
}

Int_t HMdcGeantSeg::getLastLayer12(void) const {
  if(nHitsMdc[1]) {
    for(Int_t l=5;l>=0;l--) if(geantLay[1][l]) return l+6;
  } else for(Int_t l=5;l>=0;l--) if(geantLay[0][l]) return l;
  return -1;
}

HGeantMdc* HMdcGeantSeg::getMdcLayerHit(Int_t m,Int_t l) {
  return (nMdcOk(m) && l>=0 && l<7) ? geantLay[m&1][l]:0;
}

Int_t HMdcGeantSeg::getNGMdcs(void) const {
  if(nHitsMdc[0]>0 && nHitsMdc[1]>0)   return 2;
  if(nHitsMdc[0]==0 && nHitsMdc[1]==0) return 0;
  return 1;
}
  
void HMdcGeantSeg::print(void) {
  printf("%3i) %i sec. %i seg.   GeantMdc:%4i hits (layers ",
      segNumber,sec+1,ioseg+1,getNGMdcHits());
  if(direction>0) printf("%2i -> %2i)",getFirstLayer12()+1,getLastLayer12()+1);
  else            printf("%2i -> %2i)",getLastLayer12()+1,getFirstLayer12()+1);
  if(areWiresColl) printf("   MdcCal1:%4i hits in %2i layers",
      getNDrTimes(),getNLayers());
  printf("\n");
}

Bool_t HMdcGeantSeg::getLayerHitPos(Int_t m, Int_t l, HGeomVector& hit,Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer.
  // If no GeantMdc hit in this layer and extrFlag=kTRUE
  // function does extrapolation to this layer.
  if(!nMdcOk(m) || l<0 || l>5) return kFALSE;
  Int_t modI = m & 1;
  if(geantLay[modI][l]) {
    Float_t ax,ay,tof,momLay;
    geantLay[modI][l]->getHit(ax,ay,tof,momLay);
    hit.setXYZ(ax,ay,0.);
    return kTRUE;
  }
  if(!extrFlag) return kFALSE;
  //================= If no GeantMdc hit in layer "l": =================
  if(nHitsMdc[modI]<1 || (nHitsMdc[modI]==1 && geantLay[modI][6]==0)) return kFALSE; // must be two hits at list.
  Int_t layLFirst,layLSecond,layRFirst,layRSecond;
  layLFirst=layLSecond=layRFirst=layRSecond=-1;
  for(Int_t lay=l-1;lay>=0; lay--) {
    if(l>2 && lay<3 && geantLay[modI][6]) {
      if(layLFirst<0) layLFirst=6;
      else if(layLSecond<0 && layLFirst!=6) layLSecond=6;
    }
    if(geantLay[modI][lay]) {
      if(layLFirst<0) layLFirst=lay;
      else if(layLSecond<0) layLSecond=lay;
    }
  }
  for(Int_t lay=l+1;lay<6; lay++) {
    if(l<3 && lay>2 && geantLay[modI][6]) {
      if(layRFirst<0) layRFirst=6;
      else if(layRSecond<0 && layRFirst!=6) layRSecond=6;
    }
    if(geantLay[modI][lay]) {
      if(layRFirst<0) layRFirst=lay;
      else if(layRSecond<0) layRSecond=lay;
    }
  }
  if(layLFirst<0)      calcMdcHitPos(modI,layRFirst,  layRSecond, hit,l);
  else if(layRFirst<0) calcMdcHitPos(modI,layLSecond, layLFirst,  hit,l);
  else                 calcMdcHitPos(modI,layLFirst,  layRFirst,  hit,l);
  return kTRUE;
}

Int_t HMdcGeantSeg::getFirstGeantMdcLayer(Int_t m) const {
  // Return first layer with geant hit in mdc "m".
  if(!nMdcOk(m)) return -1;
  Int_t modI = m & 1;
  for(Int_t lay=0;lay<6;lay++) if(geantLay[modI][lay]) return lay;
  return -1;
}

Int_t HMdcGeantSeg::getLastGeantMdcLayer(Int_t m) const {
  // Return first layer with geant hit in mdc "m".
  if(!nMdcOk(m)) return -1;
  Int_t modI = m & 1;
  for(Int_t lay=5;lay>=0;lay--) if(geantLay[modI][lay]) return lay;
  return -1;
}

Bool_t HMdcGeantSeg::getFirstAndLastGMdcLayers(Int_t m,Int_t& lFisrt,Int_t& lLast) const {
  // return kFALSE if "m" is not valid or one mdc geant hit exist only.
  if(!nMdcOk(m)) return kFALSE;
  Int_t modI = m & 1;
  if(nHitsMdc[modI]<2) return kFALSE;
  for(lFisrt=0;lFisrt<6;lFisrt++) if(geantLay[modI][lFisrt]) break;
  for(lLast=5; lLast>=0;lLast--)  if(geantLay[modI][lLast])  break;
  return kTRUE;
}

Bool_t HMdcGeantSeg::getMdcHitPos(Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mdc mid-plane.
  if(!nMdcOk(m)) return kFALSE;
  Int_t modI = m & 1;
  if(geantLay[modI][6]) {
    Float_t ax,ay,tof,momLay;
    geantLay[modI][6]->getHit(ax,ay,tof,momLay);
    hit.setXYZ(ax,ay,0.);
    return kTRUE;
  }
  //=============== If no GeantMdc hit in mdc mid-plane: ===============
  if(nHitsMdc[modI]<2) return kFALSE;
  Int_t layLFirst,layLSecond,layRFirst,layRSecond;
  layLFirst=layLSecond=layRFirst=layRSecond=-1;
  for(Int_t lay=2;lay>=0; lay--) if(geantLay[modI][lay]) {
    if(layLFirst<0)       layLFirst  = lay;
    else if(layLSecond<0) layLSecond = lay;
  }
  for(Int_t lay=3;lay<6; lay++) if(geantLay[modI][lay]) {
    if(layRFirst<0)       layRFirst  = lay;
    else if(layRSecond<0) layRSecond = lay;
  }
  if(layLFirst<0)      calcMdcHitPos(modI,layRFirst,  layRSecond, hit);
  else if(layRFirst<0) calcMdcHitPos(modI,layLSecond, layLFirst,  hit);
  else                 calcMdcHitPos(modI,layLFirst,  layRFirst,  hit);
  return kTRUE;
}

void HMdcGeantSeg::calcMdcHitPos(Int_t modI,Int_t lay1, Int_t lay2,HGeomVector& hit, Int_t lay) {
  // Linear extrapolation of two geant mdc hits to the mdc mid-plane.
  // Not precise solution but well enough!
  Float_t x1,y1,x2,y2,tof,momLay;
  geantLay[modI][lay1]->getHit(x1,y1,tof,momLay);
  geantLay[modI][lay2]->getHit(x2,y2,tof,momLay);
  Float_t l0 = (lay==6)  ? 2.5:lay;
  Float_t l1 = (lay1==6) ? 2.5:lay1;
  Float_t l2 = (lay2==6) ? 2.5:lay2;
  Double_t norm=(l0-(l1))/(l2-l1);
  hit.setXYZ(x1+(x2-x1)*norm,y1+(y2-y1)*norm,0.);
}

Bool_t HMdcGeantSeg::getMdcHitPosSec(Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mid-plane mdc in sector coor.sys.
  HMdcSizesCells* pMdcSizesCells=HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return kFALSE;
  if(!getMdcHitPos(m,hit)) return kFALSE;
  HMdcSizesCellsSec& pSCSec = (*pMdcSizesCells)[sec];
  HMdcSizesCellsMod& pSCMod = pSCSec[m];
  hit=pSCMod.getSecTrans()->transFrom(hit);
  return kTRUE;
}

Bool_t HMdcGeantSeg::getMdcHitPosLab(Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mid-plane mdc in lab.coor.sys..
  HMdcSizesCells* pMdcSizesCells=HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return kFALSE;
  if(!getMdcHitPos(m,hit)) return kFALSE;
  HMdcSizesCellsSec& pSCSec = (*pMdcSizesCells)[sec];
  HMdcSizesCellsMod& pSCMod = pSCSec[m];
  hit = pSCMod.getSecTrans()->transFrom(hit);
  hit = pSCSec.getLabTrans()->transFrom(hit);
  return kTRUE;
}

Bool_t HMdcGeantSeg::getLayerHitPosSec(Int_t m, Int_t l,HGeomVector& hit,Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer in sector coor.sys..
  HMdcSizesCells* pMdcSizesCells=HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return kFALSE;
  if(!getLayerHitPos(m,l,hit,extrFlag)) return kFALSE;
  HMdcSizesCellsSec&   pSCSec = (*pMdcSizesCells)[sec];
  HMdcSizesCellsLayer& pSCLay = pSCSec[m][l];
  hit = pSCLay.getSecTransGeant().transFrom(hit);
  return kTRUE;
}

Bool_t HMdcGeantSeg::getLayerHitPosLab(Int_t m, Int_t l, HGeomVector& hit,
    Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer in lab.coor.sys.
  // If no GeantMdc hit in this layer and extrFlag=kTRUE
  // function does extrapolation to this layer.
  HMdcSizesCells* pMdcSizesCells = HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return kFALSE;
  if(!getLayerHitPos(m,l,hit,extrFlag)) return kFALSE;
  HMdcSizesCellsSec&   pSCSec = (*pMdcSizesCells)[sec];
  HMdcSizesCellsLayer& pSCLay = pSCSec[m][l];
  hit = pSCLay.getSecTransGeant().transFrom(hit);
  hit = pSCSec.getLabTrans()->transFrom(hit);
  return kFALSE;
}

void HMdcGeantSeg::setStatusFlags(UChar_t& trackStatus) {
  segmentStatus = trackStatus;
  if( !areWiresColl ) return;
  if(getNLayers()<3 || getNDrTimes()<5) {
    trackStatus   &= NOT_RECONSTRUCTABLE;
    segmentStatus  = trackStatus;
    segmentStatus &= SEG_NOT_RECONS_ABLE;
  }
  if(HMdcGeantTrack::isMdcActive(sec,ioseg*2)   && (getNLayersMod(0)<3 ||
     getNDrTimesMod(0)<5)) segmentStatus &= HIT_NOT_RECONS_ABLE;
  if(HMdcGeantTrack::isMdcActive(sec,ioseg*2+1) && (getNLayersMod(1)<3 ||
     getNDrTimesMod(1)<5)) segmentStatus &= HIT_NOT_RECONS_ABLE;
}

void HMdcGeantSeg::analyseClust(HMdcClusSim* pClus,Int_t trInd,Int_t modICl) {
  // Analyse one mod. of seg. if track present in one mod. only or it is mixed cl.:
  Int_t                                         mInd = -2;
  if(pClus->getNLayersInTrack(trInd,0)==0)      mInd = 1;      //Track pres.in one mdc
  else if(pClus->getNLayersInTrack(trInd,1)==0) mInd = 0;      //Track pres.in one mdc
  else if(modICl>=0)                            mInd = modICl; //Mixed cl.finder

  UChar_t status = segmentStatus;
  if( !pClus->isSegmentAmpCut() ) {
    if(pClus->getMinCl(0) > getNLayersMod(0) ||
       pClus->getMinCl(1) > getNLayersMod(1))   status &= CLFLEVEL_TOO_HIGH;
  } else if(pClus->getMinCl() > getNLayers())   status &= CLFLEVEL_TOO_HIGH;

//  if(pClus->getNLayOrientation(trInd) < 2)      status &= FAKE_CONTRIBUTION;
  if( !pClus->is40degCross(trInd) )             status &= FAKE_CONTRIBUTION;
  else if(mInd < 0) {
         if(!isSegClusBetter(pClus,trInd))      status &= FAKE_CONTRIBUTION;
  } else if(!isModClusBetter(pClus,trInd,mInd)) status &= FAKE_CONTRIBUTION;
  pClus->setTrackStatus(trInd,status);
}

Bool_t HMdcGeantSeg::isSegClusBetter(HMdcClusSim* pClus,Int_t trInd) {
  Int_t flag = setSegClustPos(pClus,trInd);
  if(flag) {
    if(flag == 1) return kTRUE;   // Address mismatching or no HMdcSizesCells object
    if(flag == 2) return kFALSE;  // Number layers in cluster <2 from track "trInd"
    if(flag > 2)  return kTRUE;   // Not enough geant mdc hits and others
  }

  Int_t   nLayersCl = pClus->getNLayersInTrack(trInd);
  UChar_t nWiresCl  = pClus->getNTimesInTrack(trInd);
  if(pClusBest != NULL) {
    if(nLayers > nLayersCl) return kFALSE; // fake contribution
    if(nLayers == nLayersCl) {
      if(nWires > nWiresCl) return kFALSE; // fake contrirution
      if(nWires == nWiresCl) {
        Float_t dXn = pClus->dX(trInd);    // Xgeant - pClus->getX();
        Float_t dYn = pClus->dY(trInd);    // Ygeant - pClus->getY();
        if(dXn/2.7*dXn/2.7+dYn*dYn > dX/2.7*dX/2.7+dY*dY) return kFALSE;
      }
    }
    pClusBest->setTrackStatus(trIndBest,pClusBest->getTrackStatus(trIndBest) & FAKE_CONTRIBUTION);
  } else if(pClusBestCh[0] != NULL || pClusBestCh[1] != NULL) {
    // Comparing segement and chamber/mixed clusters:
    for(Int_t m=0;m<2;m++) if(pClusBestCh[m]) {
      if(nLayersCh[m] > nLayersCl)                        return kFALSE;
      if(nLayersCh[m]==nLayersCl && nWiresCh[m]>nWiresCl) return kFALSE;
    }
    // Segmet cluster better:
    if(pClusBestCh[0]) pClusBestCh[0]->setTrackStatus(trIndBestCh[0],
        pClusBestCh[0]->getTrackStatus(trIndBestCh[0]) & FAKE_CONTRIBUTION);
    if(pClusBestCh[1]) pClusBestCh[1]->setTrackStatus(trIndBestCh[1],
        pClusBestCh[1]->getTrackStatus(trIndBestCh[1]) & FAKE_CONTRIBUTION);
    pClusBestCh[0] = pClusBestCh[1] = 0;
  }

  pClusBest = pClus;
  trIndBest = trInd;
  nLayers   = nLayersCl;
  nWires    = nWiresCl;
  dX        = pClus->dX(trInd); //x - pClus->getX();
  dY        = pClus->dY(trInd); //y - pClus->getY();
  return kTRUE;
}

Bool_t HMdcGeantSeg::isModClusBetter(HMdcClusSim* pClus,Int_t trInd,Int_t mInd) {
  Int_t flag = setModClustPos(pClus,trInd,mInd);
  if(flag) {
    if(flag == 1) return kTRUE;   // Address mismatching or no HMdcSizesCells object
    if(flag == 2) return kFALSE;  // Number layers in cluster <2 from track "trInd"
    if(flag >  2) return kTRUE;   // Not enough geant mdc hits and others
  }

  Int_t   nLayersCl = pClus->getNLayersInTrack(trInd);
  UChar_t nWiresCl  = pClus->getNTimesInTrack(trInd);

  if(pClusBest) { // if exist track segment clusters:
    if(nLayers  > nLayersCl)                    return kFALSE;
    if(nLayers==nLayersCl && nWires>=nWiresCl)  return kFALSE;
    pClusBest->setTrackStatus(trIndBest,pClusBest->getTrackStatus(trIndBest) & FAKE_CONTRIBUTION);
    pClusBest = 0;
  } else if(pClusBestCh[mInd]) {
    if(nLayersCh[mInd] > nLayersCl) return kFALSE; // fake contribution
    if(nLayersCh[mInd] == nLayersCl) {
      if(nWiresCh[mInd] > nWiresCl) return kFALSE; // fake contrirution
      if(nWiresCh[mInd] == nWiresCl) {
        Float_t dXn = pClus->dX(trInd); //x - pClus->getX();
        Float_t dYn = pClus->dY(trInd); //y - pClus->getY();
        if(dXn/2.7*dXn/2.7+dYn*dYn > 
           dXCh[mInd]/2.7*dXCh[mInd]/2.7+dYCh[mInd]*dYCh[mInd]) return kFALSE;
      }
    }
    pClusBestCh[mInd]->setTrackStatus(trIndBestCh[mInd],
        pClusBestCh[mInd]->getTrackStatus(trIndBestCh[mInd])&FAKE_CONTRIBUTION);
  }
  pClusBestCh[mInd] = pClus;
  trIndBestCh[mInd] = trInd;
  nLayersCh[mInd]   = nLayersCl;
  nWiresCh[mInd]    = nWiresCl;
  dXCh[mInd]        = pClus->dX(trInd); //x - pClus->getX();
  dYCh[mInd]        = pClus->dY(trInd); //y - pClus->getY();

  if(pClusBestCh[0] && pClusBestCh[1]) {
     // Test for the same wires in two clusters:
     HMdcList12GroupCells lId1,lId2;
     compare(pClusBestCh[0],0,11, &lId1);
     compare(pClusBestCh[1],0,11, &lId2);
     Int_t n1  = lId1.getNDrTimes();
     Int_t n2  = lId2.getNDrTimes();
     Int_t n12 = lId1.nIdentDrTimes(&lId2);
     if(n1==n12 || n2==n12) {
       Int_t mRem = 0;
       if(n1 == n2) {
         if(dXCh[1]/2.7*dXCh[1]/2.7+dYCh[1]*dYCh[1] >
            dXCh[0]/2.7*dXCh[0]/2.7+dYCh[0]*dYCh[0]) mRem = 1;
       } else if(n2==n12) mRem = 1;
       if(mRem == mInd) {
         pClusBestCh[mRem]=0;
         return kFALSE;
       }
       pClusBestCh[mRem]->setTrackStatus(trIndBestCh[mRem],pClusBestCh[mRem]->
           getTrackStatus(trIndBestCh[mRem])&FAKE_CONTRIBUTION);
       pClusBestCh[mRem]=0;
    }
  }
  return kTRUE;
}

Int_t HMdcGeantSeg::setModClustPos(HMdcClusSim* pMdcClusSim, Int_t indtr,
    Int_t mInd) {
  // calculate geant track position on the project plane for chamber claster
  // return x,y and flag
  // flag = 0 - ok.
  // flag = 1 - address mismatching or no HMdcSizesCells object
  // flag = 2 - number layers in cluster <2 from track "trInd"
  // flag = 3 - not enough geant mdc hits and others
  if(pMdcClusSim->getSec() != sec) return 1;
  if(pMdcClusSim->getIOSeg() != ioseg) return 1;
  HMdcSizesCells* pMdcSizesCells=HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return 1;
  if(pMdcClusSim->getNLayersInTrack(indtr,mInd) < 2) return 2; // ??????????????????

  HGeomVector hit;
  if(!getMdcHitPosSec(mInd+ioseg*2,hit)) return 3;
  Float_t x,y;
  pMdcClusSim->calcIntersection(hit,x,y);
  pMdcClusSim->setXYGeant(indtr,x,y);
  return 0;
// //----------------END???????????
// 
//   HGeomVector hit1;
//   HGeomVector hit2;
// Float_t dX1,dY1,dX2,dY2,dX4,dY4,dX,dY;
// dX1=dY1=dX2=dY2=dX=dY=-100000.;
// 
// 
//   if(mod>=0) {
//     if(!getMdcHitPosSec(mod,hit1)) return 3;
//     pMdcClusSim->calcIntersection(hit1,x,y);
// pMdcClusSim->print();
// printf("+++ dX3=%.1f dY3=%.1f \n",pMdcClusSim->getX()-x,pMdcClusSim->getY()-y);
//     return 0;
//   }
// 
// 
// pMdcClusSim->print();
// printf("+++=0! nLaysM1=%i nLaysM2=%i\n",pMdcClusSim->getNLayersInTrack(indtr,0),pMdcClusSim->getNLayersInTrack(indtr,1));
// //?    Int_t m = (nLaysM1>0) ? ioseg*2 : ioseg*2+1;
//     Int_t m = mInd + ioseg*2;
//     Int_t l1,l2;
//     if(!getFirstAndLastGMdcLayers(m,l1,l2)) return 3;
// printf("+++ m=%i l1=%i l2=%i\n",m,l1,l2);
//     if(!getLayerHitPosSec(m,l1,hit1,kFALSE)) return 3;
//     if(!getLayerHitPosSec(m,l2,hit2,kFALSE)) return 3;
// Float_t x1,y1,x2,y2;
// pMdcClusSim->calcIntersection(hit1,x1,y1);
// pMdcClusSim->calcIntersection(hit2,x2,y2);
// dX1=pMdcClusSim->getX()-x1;
// dX2=pMdcClusSim->getX()-x2;
// dY1=pMdcClusSim->getY()-y1;
// dY2=pMdcClusSim->getY()-y2;
// 
// //? if(nLaysM2==0) {if(!getMdcHitPosSec(ioseg*2,hit)) return 3;}
// //? else if(!getMdcHitPosSec(ioseg*2+1,hit)) return 3;
// if(!getMdcHitPosSec(m,hit)) return 3;
// pMdcClusSim->calcIntersection(hit,x1,y1);
// dX4=pMdcClusSim->getX()-x1;
// dY4=pMdcClusSim->getY()-y1;
// printf("+++ dX4=%.1f dY4=%.1f \n",pMdcClusSim->getX()-x1,pMdcClusSim->getY()-y1);
// 
//   pMdcClusSim->calcIntersection(hit1,hit2,x,y);
// dX=pMdcClusSim->getX()-x;
// dY=pMdcClusSim->getY()-y;
// printf("***  dX=%.1f  dY=%.1f\n",dX,dY);
// if(dX1>-90000.) printf("+++ dX1=%.1f dY1=%.1f \n+++ dX2=%.1f dY2=%.1f\n",dX1,dY1,dX2,dY2);
//   return 0;
}

Int_t HMdcGeantSeg::setSegClustPos(HMdcClusSim* pMdcClusSim, Int_t indtr) {
  // calculate geant track position on the project plane for segment claster
  // return flag:
  // flag = 0 - ok.
  // flag = 1 - address mismatching or no HMdcSizesCells object
  // flag = 2 - number layers in track "trInd" <2
  // flag = 3 - not enough geant mdc hits and others
  if(pMdcClusSim->getSec() != sec) return 1;
  if(pMdcClusSim->getIOSeg() != ioseg) return 1;
  HMdcSizesCells* pMdcSizesCells=HMdcSizesCells::getExObject();
  if(pMdcSizesCells==0) return 1;

  if(pMdcClusSim->getNLayersInTrack(indtr) < 2) return 2; // ??????????????????

  HGeomVector hit1;
  HGeomVector hit2;
  if(!getMdcHitPosSec(ioseg*2,  hit1)) {
    if(nHitsMdc[0] <= 0) return 3;
    for(Int_t l=5;l>=0;l--) if(geantLay[0][l]) {
      if(!getLayerHitPosSec(ioseg*2,l,hit1,kFALSE)) return 3;
      break;
    }
  }
  if(!getMdcHitPosSec(ioseg*2+1,hit2)) {
    if(nHitsMdc[1] <= 0) return 3;
    for(Int_t l=0;l<6;l++) if(geantLay[1][l]) {
      if(!getLayerHitPosSec(ioseg*2+1,l,hit2,kFALSE)) return 3;
      break;
    }
  }
  Float_t x,y;
  pMdcClusSim->calcIntersection(hit1,hit2,x,y);
  pMdcClusSim->setXYGeant(indtr,x,y);

// pMdcClusSim->print();
// Float_t dX1,dY1,dX2,dY2,dX,dY;
// dX1=dY1=dX2=dY2=dX=dY=-100000.;
// Float_t x1,y1,x2,y2;
// pMdcClusSim->calcIntersection(hit1,x1,y1);
// pMdcClusSim->calcIntersection(hit2,x2,y2);
// dX1=pMdcClusSim->getX()-x1;
// dX2=pMdcClusSim->getX()-x2;
// dY1=pMdcClusSim->getY()-y1;
// dY2=pMdcClusSim->getY()-y2;
// 
// dX=pMdcClusSim->getX()-x;
// dY=pMdcClusSim->getY()-y;
// printf("***  dX=%.1f  dY=%.1f\n",dX,dY);
// if(dX1>-90000.) printf("+++ dX1=%.1f dY1=%.1f \n+++ dX2=%.1f dY2=%.1f\n",dX1,dY1,dX2,dY2);

  return 0;
}

void HMdcGeantSeg::analyseSeg(HMdcSegSim* pSeg,Int_t trInd,HMdcList12GroupCells& wrList,
                              Float_t ch2,Int_t nR) {
  UChar_t status = segmentStatus | SEGS_MATH_OK; // No seg.match. info in MdcSegSim any more!
  if( !HMdcBArray::is40DegWireCross(wrList.getListLayers()) )    status &= FAKE_CONTRIBUTION;
  else if( nR >=  pSeg->getNTimes(trInd) )                       status &= FAKE_CONTRIBUTION;
  else if( !isSegBetter(pSeg,trInd,wrList.getNLayers(),ch2,nR) ) status &= FAKE_CONTRIBUTION;
  pSeg->setTrackStatus(trInd,status);
}

Bool_t HMdcGeantSeg::isSegBetter(HMdcSegSim* pSeg,Int_t trInd,Int_t nLayersN,Float_t chi2N,Int_t nR) {
  if(pSeg->getChi2() < 0.) chi2N = -1.;
  UChar_t nWiresN  = pSeg->getNTimes(trInd);
 
  if(pSegBest != NULL) {
    if(chi2N < 0. && chi2 >= 0.)   return kFALSE;   // fake contribution
    else   if(chi2N <  0. && chi2 <  0.) {
      if(sNLayers > nLayersN) return kFALSE;   // fake contribution
      if(sNLayers == nLayersN) {
        if(sNWires-nWiresN > 1) return kFALSE; // fake contrirution
        if(sNWires-nWiresN >= -1 && nR > nRestWires) return kFALSE;
      }
    } else if(chi2N >= 0. && chi2 >= 0.) {
      if(sNLayers > nLayersN) return kFALSE;   // fake contribution
      if(sNLayers == nLayersN) {
        if(sNWires-nWiresN > 1) return kFALSE; // fake contrirution
        if(sNWires-nWiresN >= -1) {
          if(nR-nRestWires >  1) return kFALSE;
          if(nR-nRestWires > -1 && chi2N > chi2) return kFALSE;
        }
      }
    }
    pSegBest->setTrackStatus(sTrIndBest,pSegBest->getTrackStatus(sTrIndBest) & FAKE_CONTRIBUTION);
  }
  pSegBest   = pSeg;
  sTrIndBest = trInd;
  sNLayers   = nLayersN;
  sNWires    = nWiresN;
  chi2       = chi2N;
  nRestWires = nR;
  return kTRUE;
}

//===========================================================================

Bool_t  HMdcGeantTrack::mdcSetup[6][4]=
    {{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1},{1,1,1,1}};

HMdcGeantTrack::HMdcGeantTrack(void) {
  // constructor for the case when geantSegment array created in HMdcGeantTrack
  arrSizeLoc    = 0;
  arrSize       = &arrSizeLoc;
  arrOffset     = 0;
  arrGlobOffset = 0;
  iterMdcCal1   = 0;
  segments      = new TObjArray();
  clear();
}

HMdcGeantTrack::HMdcGeantTrack(TObjArray* arr,Int_t* size,Int_t* offst) {
  // constructor for the case when geantSegment array created out of HMdcGeantTrack
  // parameters:
  // arr - pointer to TObjArray for HMdcGeantSeg objects
  // size - pointer to Int_t value which keep num. of not empty slots
  //        in TObjArray (should not has empty slots between not empty slots)
  // offset - pointer to Int_t value which keep num. of HMdcGeantSeg objects
  //          already used in TObjArray for previous tracks
  segments      = arr;
  arrSize       = size;
  arrOffset     = 0;
  arrGlobOffset = offst;
  pMdcCal1Cat   = NULL;
  iterMdcCal1   = NULL;
  clear();
}

HMdcGeantTrack::~HMdcGeantTrack(void) {
  if(iterMdcCal1) {
    delete iterMdcCal1;
    iterMdcCal1 = 0;
  }
  if(segments && arrGlobOffset==0) {
    segments->Delete();
    delete segments;
    segments=0;
  }
}

void HMdcGeantTrack::clear(void) {
  pGeantKine     = 0;
  debugPrintFlag = kFALSE;
  setDefault();
}

void HMdcGeantTrack::setDefault(void) {
  nGeantMdcHits = 0;
  nSegments     = 0;
  nWSegments    = 0;
  mdcWSector    = -1;
  mdcSector     = -1;
  nIOSeg        = 0;
  nWIOSeg       = 0;
  nWires        = -1;
  directionFlag = kTRUE;
  geantBugFlag  = kFALSE;
  trackStatus   = DEFAULT_SET;
  userFlag      = 0;
}

void HMdcGeantTrack::testMdcSetup(void) {
  if(gHades==0) return;
  for(Int_t s=0;s<6;s++) for(Int_t m=0;m<4;m++) 
      mdcSetup[s][m] = HMdcGetContainers::getObject()->isModActive(s,m);
}

Bool_t HMdcGeantTrack::setMdcCal1Cat(void) {
  HCategory* cat = HMdcGetContainers::getObject()->getCatMdcCal1();
  if(pMdcCal1Cat && pMdcCal1Cat==cat) return kTRUE;
  pMdcCal1Cat=cat;
  if(pMdcCal1Cat==0) return kFALSE;
  if(iterMdcCal1) delete iterMdcCal1;
  iterMdcCal1 = (HIterator*)pMdcCal1Cat->MakeIterator();
  return kTRUE;
}

Short_t HMdcGeantTrack::testGeantTrack(Int_t trNum) {
  // return number of mdc segments
  HCategory* pGeantKineCat = HMdcGetContainers::getObject()->getCatGeantKine();
  if(pGeantKineCat) return testGeantTrack((HGeantKine*)pGeantKineCat->getObject(trNum-1));
  Error("testGeantTrack"," Can't get catGeantKine category!");
  setDefault();
  return 0;
}

Short_t HMdcGeantTrack::testGeantTrack(HGeantKine* pGK) {
  // return number of mdc segments
  pGeantKine = pGK;
  if(pGeantKine==0) {
    setDefault();
    Error("testGeantTrack"," pointer of HGeantKine object = 0!");
    return 0;
  }
  testHitsInDetectors();
  return testMdcHits();
}

void HMdcGeantTrack::testHitsInDetectors(void) {
  trackNumber  = pGeantKine->getTrack();
  mom          = pGeantKine->getTotalMomentum();
  isInMdcFlag  = pGeantKine->getFirstMdcHit()    >= 0;
  isInRichFlag = pGeantKine->getFirstRichHit()   >= 0;
  isInMetaFlag = pGeantKine->getFirstTofHit()    >= 0 ||
                 pGeantKine->getFirstShowerHit() >= 0 ||
                 pGeantKine->getFirstRpcHit()    >= 0;
  if(!isInMetaFlag) trackStatus &= NO_META_HIT;
  if(!isInMdcFlag)  trackStatus &= NOT_RECONSTRUCTABLE;
}

Short_t HMdcGeantTrack::testMdcHits(void) {
  // return number of mdc segments
  setDefault();
  if(!isInMdcFlag) return 0;
  if(arrGlobOffset) arrOffset = *arrGlobOffset;

  pGeantMdcPrev     = 0;
  Char_t  modPrev   = -1;
  Char_t  lay14Prev = -1;
  Float_t tofPrev   = -1.;
  Float_t momPrev   = mom;
  geantBugFlag      = kFALSE;

  if(!addSegAtAndExpand(0)) return 0;
  HCategory* pGeantMdcCat = HMdcGetContainers::getObject()->getCatGeantMdc();
  if(pGeantMdcCat==0) {
    Error("testMdcHits","Can't get catMdcGeantRaw category!");
    return 0;
  }
  pGeantKine->setMdcCategory(pGeantMdcCat);
  pGeantKine->resetMdcIter();
  while((pGeantMdc=(HGeantMdc*)pGeantKine->nextMdcHit()) != NULL) {
    Char_t sec = pGeantMdc->getSector();
    Char_t mod = pGeantMdc->getModule();
    if(!mdcSetup[(Int_t)sec][(Int_t)mod]) continue;
    Char_t lay    = pGeantMdc->getLayer();
    Char_t lay14  = calcLay14(mod,lay);
    Char_t hitDir = HMdcGeantSeg::dirHit(pGeantMdc);
    Float_t ax,ay,tof,momLay;
    pGeantMdc->getHit(ax,ay,tof,momLay);
    Char_t  dMod  = mod-modPrev;
    Char_t  dLay  = lay14-lay14Prev;
    Float_t dMom  = momPrev-momLay;
    Bool_t isDirOk = kTRUE;
    if(modPrev>=0 && (hitDir!=dirDLayer(dLay) || hitDir!=segment->getDirection())) isDirOk = kFALSE;

    // Finding bugs in HGeantMdc:
//     if(isGeantBug1(momLay)) continue;
//     if(isGeantBug2(sec,dMod,dLay,hitDir,tof-tofPrev)) continue;
//     mayBeGeantBug(dMom);
    if( !isGeantBug1(momLay) && !isGeantBug2(sec,dMod,dLay,hitDir,tof-tofPrev)) mayBeGeantBug(dMom);

    if(debugPrintFlag && nGeantMdcHits==0) printf("---- New track: ----\n");
    if(modPrev<0 || !isDirOk || !segment->addHit(pGeantMdc)) {
      // New segment ---------------------------
      if(modPrev>=0 && !addSegment()) break;
      segment->addFirstHit(pGeantMdc,mdcSetup[(Int_t)sec],nSegments);
      segment->setTrackNumber(trackNumber);
    }
    if(lay!=6) collectWires(sec,mod,lay,tof);
    if(debugPrintFlag) printDebug(dMom,dMod);
    nGeantMdcHits++;
    pGeantMdcPrev = pGeantMdc;
    modPrev       = mod;
    lay14Prev     = lay14;
    tofPrev       = tof;
    momPrev       = momLay;
  }
  if(nGeantMdcHits>0) {
    if(segment->getNGMdcHits()>0) nSegments++;
    else nGeantMdcHits--;  // remove one mid-planeHit segment
  }
  if(nGeantMdcHits==0) isInMdcFlag=kFALSE;
  calcNSectors();
  setStatusFlags();
  return nSegments;
}

void HMdcGeantTrack::setStatusFlags(void) {
  Bool_t* mdcSp = mdcSector>=0 ? mdcSetup[(Int_t)mdcSector] : 0;
  if(geantBugFlag)                          trackStatus &= GEANT_BUG;
  if(nSegments<1 || nSegments>2)            trackStatus &= NOT_RECONSTRUCTABLE;
  else if(mdcSector < 0)                    trackStatus &= NOT_RECONSTRUCTABLE;
  else if(!directionFlag)                   trackStatus &= NOT_RECONSTRUCTABLE;
  else if(nSegments == 1) {
    if(nIOSeg==1 && (mdcSp[2] || mdcSp[3])) trackStatus &= NOT_RECONSTRUCTABLE;
    if(nIOSeg==2 && (mdcSp[0] || mdcSp[1])) trackStatus &= NOT_RECONSTRUCTABLE;
  } else if(nWires >= 0) {
    // wires are collected
    if(nWSegments<1 || nWSegments>2)        trackStatus &= NOT_RECONSTRUCTABLE;
    else if(nWires < 5*nWSegments)          trackStatus &= NOT_RECONSTRUCTABLE;
    else if(nWSegments == 1) {
      if(nWIOSeg==1) {
        if(mdcSp[2] || mdcSp[3])            trackStatus &= NOT_RECONSTRUCTABLE;
      } else if(nWIOSeg==2) {
        if(mdcSp[0] || mdcSp[1])            trackStatus &= NOT_RECONSTRUCTABLE;
      }
    }
  }
  for(Short_t ns=0;ns<nSegments;ns++) {
    HMdcGeantSeg* pGeantSeg = getSegment(ns);
    pGeantSeg->setStatusFlags(trackStatus);
  }
}

Bool_t HMdcGeantTrack::addSegAtAndExpand(Short_t segNum) {
  if(*arrSize <= segNum+arrOffset) {
    segment = new HMdcGeantSeg(segNum);
    segments->AddAtAndExpand(segment,segNum+arrOffset);
    (*arrSize)++;
  } else {
    segment = (HMdcGeantSeg*)segments->At(segNum+arrOffset);
    if(segment==0) return kFALSE;
    segment->clear(segNum);
  }
  return kTRUE;
}

Bool_t HMdcGeantTrack::addSegment(void) {
  if(segment->getNGMdcHits()>0) {
    nSegments++;
    if(!addSegAtAndExpand(nSegments)) return kFALSE;
    if(debugPrintFlag) printf("- New segment:\n");
  } else {
    if(debugPrintFlag) printDebug(pGeantMdcPrev,-1," ! mid-plane segment");
    nGeantMdcHits--;  // remove segment from one mdc mid-plane hit
  }
  return kTRUE;
}

Int_t HMdcGeantTrack::calcLay14(Int_t m, Int_t l) {
  // recalculate module and layer number in layer number from 0 to 13
  // were numb. 3 and 10 correspond mid-planes of mdc's
  Int_t ls = (m&1)*7;
  if(l<3) return l+ls;
  if(l<6) return l+1+ls;
  return 3+ls;             // mid-plane MDC
}

Bool_t HMdcGeantTrack::isGeantBug1(Float_t momLay) {
  // Finding bug in HGeantMdc. Hit momentum > vertex momentum.
  if(mom < momLay-0.001) {
    if(debugPrintFlag) printDebug(pGeantMdc,-1,"! P>Pver");
    geantBugFlag = kTRUE;
    return kTRUE;
  }
  return kFALSE;
}

Bool_t HMdcGeantTrack::isGeantBug2(Char_t sec,Char_t dMod, Char_t dLay,
    Char_t hitDir, Float_t dTof) {
  // Finding bug in HGeantMdc. Too short t.o.f. between two hits.
  if(dTof>0.006 || segment->getNGMdcHits()<=0) return kFALSE;
  if(sec==segment->getSec() && dLay==0 && dMod==0) {
//    if(hitDir != segment->getDirection() || !isGeantBug3()) return kFALSE;
    if(hitDir != segment->getDirection() /*|| !isGeantBug3()*/) return kFALSE;
  } else if(abs(dLay) <= 1) return kFALSE;
  if(debugPrintFlag) printDebug(pGeantMdc,-1,"! dTof<0.006");
  geantBugFlag = kTRUE;
  return kTRUE;
}

Bool_t HMdcGeantTrack::mayBeGeantBug(Float_t dMom) {
  // Finding bug in HGeantMdc. Hit momentum > prev.hit momentum.
  if(dMom > -0.01) return kFALSE;
  geantBugFlag = kTRUE;
  return kTRUE;
}

Bool_t HMdcGeantTrack::isGeantBug3(void) {
  // removing fake HGeantMdc hits from list
  Float_t x1,y1,tof1,p1;
  pGeantMdcPrev->getHit(x1,y1,tof1,p1);
  Float_t x2,y2,tof2,p2;
  pGeantMdc->getHit(x2,y2,tof2,p2);
  if(fabs(x1-x2) > 4.) return kFALSE;
  if(fabs(y1-y2) > 4.) return kFALSE;
  Char_t lay=pGeantMdc->getLayer();
  if(lay!=6) {
    Int_t nwires1 = segment->getNDrTimes();
    Int_t nwires2 = collectWires(pGeantMdc->getSector(),pGeantMdc->getModule(),lay,tof2);
    if(nwires2 > nwires1) {
      pGeantMdcPrev = pGeantMdc;
      segment->setAnotherHit(pGeantMdc);
    }
  }
  return kTRUE;
}

Int_t HMdcGeantTrack::collectWires(Char_t sec, Char_t mod, Char_t lay, Float_t atof) {
  Int_t nwires = 0;
  HMdcEvntListCells* pMdcListCells = HMdcEvntListCells::getExObject();
  if(pMdcListCells) {
    HMdcLayListCells& pLayListCells=(*pMdcListCells)[sec][mod][lay];
    Int_t cell = -1;
    Int_t track1,track2;
    Float_t tof1,tof2;
    UChar_t times;
    while( (times=pLayListCells.nextCellSim(cell,track1,track2,tof1,tof2)) ) {
      if(track1!=trackNumber || fabs(atof-tof1)>0.0005) times &= 2;
      if(track2!=trackNumber || fabs(atof-tof2)>0.0005) times &= 1;
      if(times==0) continue;
      //Int_t notAdded = segment->setTime(lay+(mod&1)*6,cell,times);
      segment->setTime(lay+(mod&1)*6,cell,times);
// Int_t gtr=pLayListCells.getGeantTrack(cell,1);
// if(notAdded!=0) {
// //  HGeantKine *kine = (HGeantKine*)pGeantKineCat->getObject(gtr-1);
//   printf("nn=%i gntr=%i pid=%i %iL %ic atof=%f\n",notAdded,gtr,pGeantKine->getID(),lay+(mod&1)*6,cell,atof);
// }
      nwires++;
    }
    segment->setWiresAreColl();
//if(nwires>0) printf("************** nwires=%i\n",nwires);
    return nwires;
  } else if(setMdcCal1Cat()) {
    locMdcCal1.set(3,sec,mod,lay);
    iterMdcCal1->Reset();
    iterMdcCal1->gotoLocation(locMdcCal1);
    HMdcCal1Sim* pMdcCal1Sim;
    while((pMdcCal1Sim=(HMdcCal1Sim*)iterMdcCal1->Next())) {
      Int_t nHits=pMdcCal1Sim->getNHits();
      if(nHits==0) continue;
      UChar_t times=(nHits>0) ? 1:-nHits;
      if(times==2) times=3;
      if((times&1) && (pMdcCal1Sim->getStatus1()<=0 ||
          pMdcCal1Sim->getNTrack1()!=trackNumber ||
          fabs(atof-pMdcCal1Sim->getTof1())>0.0005)) times&=2;
      if((times&2) && (pMdcCal1Sim->getStatus2()<=0 ||
          pMdcCal1Sim->getNTrack2()!=trackNumber ||
          fabs(atof-pMdcCal1Sim->getTof2())>0.0005)) times&=1;
      if(times==0) continue;
      nwires++;
      segment->setTime(lay+(mod&1)*6,pMdcCal1Sim->getCell(),times);
    }
    segment->setWiresAreColl();
  } else segment->setWiresAreNotColl();
  return nwires;
}

void HMdcGeantTrack::unsetMdc(Int_t s, Int_t m) {
  if(s>=0 && s<6 && m>=0 && m<4) mdcSetup[s][m]=kFALSE;
}

void HMdcGeantTrack::print(void) {
  // nseg<0 - print all segments
  Float_t ax,ay,az;
  pGeantKine->getVertex(ax,ay,az);
  printf("Tr.%i  %s p=%.3f MeV. Vertex=%.1f/%.1f/%.1f ",trackNumber,
      HPhysicsConstants::pid(pGeantKine->getID()),
      pGeantKine->getTotalMomentum(),ax,ay,az);
  if(pGeantKine->isPrimary()) printf("Primary. ");
  else printf("Secondary. ");
  if(isInMdcFlag||isInRichFlag||isInMetaFlag) {
    printf("Present in");
    if(isInRichFlag) printf(" RICH");
    if(isInMdcFlag)  printf(" MDC");
    if(isInMetaFlag) printf(" META");
    printf("\n");
  } else printf("No hits.\n");
  if(nSegments>0) {
    printf(" MDC                 GeantMdc:%3i segments ",nSegments);
    if(mdcSector>=0) printf("in sector %i ",mdcSector+1);
    else  printf("in %i sectors",-mdcSector-1);
    printf("     MdcCal1:");
    if(HMdcEvntListCells::getExObject() || pMdcCal1Cat) {
      if(nWSegments==0) printf(" no hits");
      else {
        printf("%3i segments ",nWSegments);
        if(mdcWSector>=0) printf("in sector %i",mdcWSector+1);
        else printf("in %i sectors",-mdcWSector-1);
      }
    } else printf(" hits not collected");
    printf("\n");
    for(Short_t ns=0;ns<nSegments;ns++)
      ((HMdcGeantSeg*)segments->At(ns+arrOffset))->print();
  }
  if(userFlag) printf("  userFlag=%i\n",userFlag);
  printf("\n");
}

void HMdcGeantTrack::printDebug(Float_t dMom,const Char_t dMod) {
  if(dMom < -0.01) printDebug(pGeantMdc,nGeantMdcHits,"! P>Pprev.");
  else if(dMom>50. && dMod==0 && segment->getNGMdcHits()>1)
      printDebug(pGeantMdc,nGeantMdcHits,"! Pprev.-P>50");
  else printDebug(pGeantMdc,nGeantMdcHits);
}

void HMdcGeantTrack::printDebug(HGeantMdc* pGeantMdc,Int_t i,const Char_t* st) {
  Int_t sec = pGeantMdc->getSector();
  Int_t mod = pGeantMdc->getModule();
  Int_t lay = pGeantMdc->getLayer();
  Int_t lay14=calcLay14(mod,lay);
  Float_t ath, aph;
  pGeantMdc->getIncidence(ath,aph);
  Float_t ax,ay,tof,momLay;
  pGeantMdc->getHit(ax,ay,tof,momLay);
  if(i>=0) printf("%3i)",i);
  else printf("Del.");
  printf(" %3itr. %is %im %2il  l14=%2i t=%7.4f th=%7.3f ph=%7.3f",
   pGeantMdc->getTrack(),sec,mod,lay,lay14,tof,ath,aph);
  printf(" dr=%2i x=%7.1f y=%7.1f p=%.3f",HMdcGeantSeg::dirTheta(ath),
      ax,ay,momLay);
  if(st) printf(" %s",st);
  printf("\n");
}

Int_t HMdcGeantTrack::getNSegsInMdc(Int_t m,Int_t sec) {
  // return number of segments where mdc "m" present in sector "sec"
  // if "sec"<0 - return number of segments where mdc "m" present
  if(!nMdcOk(m) || nSegments==0) return 0;
  Int_t segmdc = m >> 1;
  Int_t nSeg   = 0;
  for(Short_t ns=0; ns<nSegments; ns++) {
    HMdcGeantSeg * geantSeg =(HMdcGeantSeg*)segments->At(ns+arrOffset);
    if(geantSeg->getIOSeg() != segmdc) continue;
    if(sec>=0 && geantSeg->getSec() != sec) continue;
    if(geantSeg->getNGMdcHits(m)) nSeg++;
  }
  return nSeg;
}

Int_t HMdcGeantTrack::getNSegments(Int_t seg,Int_t sec) {
  // return number of segments where mdc "m" present in sector "sec"
  // if "sec"<0 - return number of segments where mdc "m" present
  if(seg<0 || seg>1 || nSegments==0) return 0;
  Int_t nSeg=0;
  for(Short_t ns=0; ns<nSegments; ns++) {
    HMdcGeantSeg * geantSeg =(HMdcGeantSeg*)segments->At(ns+arrOffset);
    if(geantSeg->getIOSeg()  != seg) continue;
    if(sec>=0 && geantSeg->getSec() != sec)    continue;
    nSeg++;
  }
  return nSeg;
}

HGeantMdc* HMdcGeantTrack::getMdcMidPlaneHit(Int_t nseg,Int_t mod) {
  if(nSegOk(nseg)) return
    ((HMdcGeantSeg*)segments->At(nseg+arrOffset))->getMdcMidPlaneHit(mod);
  return 0;
}

HGeantMdc* HMdcGeantTrack::getMdcLayerHit(Int_t nseg,Int_t mod,Int_t lay) {
  if(nSegOk(nseg)) return
    ((HMdcGeantSeg*)segments->At(nseg+arrOffset))->getMdcLayerHit(mod,lay);
  return 0;
}

Bool_t HMdcGeantTrack::getMdcHitPos(Int_t ns, Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mid-plane mdc.
  if(nSegOk(ns)) return
    ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getMdcHitPos(m,hit);
  return kFALSE;
}

Bool_t HMdcGeantTrack::getLayerHitPos(Int_t ns, Int_t m, Int_t l,
    HGeomVector& hit, Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer.
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->
      getLayerHitPos(m,l,hit,extrFlag);
  return kFALSE;
}

Bool_t HMdcGeantTrack::getMdcHitPosSec(Int_t ns, Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mid-plane mdc in sector coor.sys.
  if(nSegOk(ns)) return
      ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getMdcHitPosSec(m,hit);
  return kFALSE;
}

Bool_t HMdcGeantTrack::getMdcHitPosLab(Int_t ns, Int_t m, HGeomVector& hit) {
  // Return geant "hit" position in mid-plane mdc in lab.coor.sys..
  if(nSegOk(ns)) return
      ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getMdcHitPosLab(m,hit);
  return kFALSE;
}

Bool_t HMdcGeantTrack::getLayerHitPosSec(Int_t ns, Int_t m, Int_t l,
    HGeomVector& hit, Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer in sector coor.sys..
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->
      getLayerHitPosSec(m,l,hit,extrFlag);
  return kFALSE;
}

Bool_t HMdcGeantTrack::getLayerHitPosLab(Int_t ns, Int_t m, Int_t l,
    HGeomVector& hit, Bool_t extrFlag) {
  // Return geant "hit" position in mdc layer in lab.coor.sys.
  // If no GeantMdc hit in this layer and extrFlag=kTRUE
  // function does extrapolation to this layer.
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->
      getLayerHitPosLab(m,l,hit,extrFlag);
  return kFALSE;
}

Char_t HMdcGeantTrack::getSegDirection(Int_t ns) {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->
      getDirection();
  return 0;
}

Char_t HMdcGeantTrack::getSector(Int_t ns)       {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getSec();
  return -1;
}

Char_t HMdcGeantTrack::getIOSeg(Int_t ns)        {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getIOSeg();
  return -1;
}

Char_t HMdcGeantTrack::getFirstLayer12(Int_t ns)   {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getFirstLayer12();
  return -1;
}

Char_t HMdcGeantTrack::getLastLayer12(Int_t ns)    {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getLastLayer12();
  return -1;
}

Char_t HMdcGeantTrack::getNGMdcHits(Int_t ns)    {
  if(nSegOk(ns)) return ((HMdcGeantSeg*)segments->At(ns+arrOffset))->getNGMdcHits();
  return 0;
}

void HMdcGeantTrack::calcNSectors(void) {
  if(nSegments==0) return;
  Short_t secList = 0;
  Int_t sec       = -1;
  for(Short_t ns=0;ns<nSegments;ns++) {
    HMdcGeantSeg* geantSeg =(HMdcGeantSeg*)segments->At(ns+arrOffset);
    if(geantSeg->getDirection()<1) directionFlag=kFALSE;
    nIOSeg |= (1<<geantSeg->getIOSeg());
    sec=geantSeg->getSec();
    Short_t tst=1<<sec;
    if(secList&tst) continue;
    secList |= tst;
    mdcSector--;
  }
  if(mdcSector==-2) mdcSector=sec;

  if(HMdcEvntListCells::getExObject()==0 && pMdcCal1Cat==0) return;
  secList = 0;
  sec     = -1;
  nWires  = 0;
  for(Short_t ns=0;ns<nSegments;ns++) {
    HMdcGeantSeg* geantSeg =(HMdcGeantSeg*)segments->At(ns+arrOffset);
    Int_t nDrTm=geantSeg->getNDrTimes();
    if(nDrTm==0) continue;
    nWires += nDrTm;
    nWIOSeg |= (1<<geantSeg->getIOSeg());
    sec=geantSeg->getSec();
    Short_t tst=1<<sec;
    nWSegments++;
    if(secList&tst) continue;
    secList |= tst;
    mdcWSector--;
  }
  if(mdcWSector==-2) mdcWSector=sec;
}

HMdcGeantTrack* HMdcGeantEvent::next(Int_t& i) {
  i=(i<0) ? 0:i+1;
  if(i>=nTracks || nTracks==0) return 0;
  return (HMdcGeantTrack*)At(i);
}

HMdcGeantTrack* HMdcGeantEvent::getGeantTrack(Int_t trackNum) {
  if(nTracks <= 0 || trackNum<=0) return 0;
  Int_t nabove = nTracks+1;
  Int_t nbelow = 0;
  while(nabove-nbelow > 1) {
     Int_t middle = (nabove+nbelow) >> 1;
     HMdcGeantTrack* pTrack=(HMdcGeantTrack*)At(middle-1);
     Int_t track = pTrack->getTrack();
     if(trackNum == track) return pTrack;
     if(trackNum  < track) nabove = middle;
     else                  nbelow = middle;
  }
  return 0;
}

//===========================================================================

HMdcGeantEvent* HMdcGeantEvent::pGlobalGEvent=0;

HMdcGeantEvent* HMdcGeantEvent::getObject(Bool_t& isCreated) {
  // In task call this function not earlier then in your init function.
  // In macro first time call this function after
  // HMdcGetContainers::setCatGeantKine(HCategory* cat);
  // HMdcGetContainers::setCatGeantMdc(HCategory* cat);
  isCreated = kFALSE;
  if(pGlobalGEvent==0) {
    if(HMdcGetContainers::getObject()->getCatGeantKine() == 0) return 0;
    if(HMdcGetContainers::getObject()->getCatGeantMdc()  == 0) return 0;
    pGlobalGEvent = new HMdcGeantEvent();
    isCreated     = kTRUE;
  }
  return pGlobalGEvent;
}

void HMdcGeantEvent::deleteCont(void) {
  if(pGlobalGEvent) delete pGlobalGEvent;
  pGlobalGEvent = 0;
}

HMdcGeantEvent::HMdcGeantEvent(void) {
  size           = 0;
  nTracks        = 0;
  pGeantKineCat  = 0;
  pGeantMdcCat   = 0;
//  pMdcCal1Cat    = 0;
  pMdcListCells  = 0;
  isMdcLCellsOwn = kFALSE;
  debugPrintFlag = kFALSE;
  nGSegments     = 0;
  sizeGSegArr    = 0;
  isMdcLCellsOwn = kFALSE;
  HMdcGeantTrack::testMdcSetup();
}

HMdcGeantEvent::~HMdcGeantEvent(void) {
  // destructor
  Delete();
  HMdcEvntListCells::deleteCont();
  geantSegments.Delete();
}

Bool_t HMdcGeantEvent::setGeantKineCat(void) {
  HCategory* cat = HMdcGetContainers::getObject()->getCatGeantKine();
  if(pGeantKineCat && cat==pGeantKineCat) return kTRUE;
  pGeantKineCat=cat;
  if(pGeantKineCat==0) {
    Error("setGeantKineCat"," Can't get catGeantKine category!");
    return kFALSE;
  }
  return kTRUE;
}

Bool_t HMdcGeantEvent::collectTracks(void) {
  nTracks      = 0;
  nGSegments   = 0;
  geantBugFlag = kFALSE;
  xVertex = yVertex = zVertex = -1000.;
  targNum = -1;
  Bool_t setVertex = kTRUE;
  if(pMdcListCells == NULL) pMdcListCells = HMdcEvntListCells::getObject(isMdcLCellsOwn);
  if(isMdcLCellsOwn) pMdcListCells->collectWires();
  if(!setGeantKineCat()) return kFALSE;
  HMdcGeantTrack* pGeantTrack;
  Int_t nTrks = pGeantKineCat->getEntries();
  for(Int_t trk=0;trk<nTrks;trk++) {
    HGeantKine* pGeantKine = (HGeantKine*)pGeantKineCat->getObject(trk);
    if(pGeantKine == NULL) continue;
    if(setVertex && pGeantKine->isPrimary()) {
      pGeantKine->getVertex(xVertex,yVertex,zVertex);
      HMdcSizesCells* pMdcSizesCells = HMdcSizesCells::getExObject();
      if(pMdcSizesCells) targNum = pMdcSizesCells->calcTargNum(zVertex);
      setVertex = kFALSE;
    }
    if(nTracks<size) pGeantTrack=(HMdcGeantTrack*)At(nTracks);
    else {
      pGeantTrack = new HMdcGeantTrack(&geantSegments,&sizeGSegArr,&nGSegments);
      AddAtAndExpand(pGeantTrack, nTracks);
      size++;
    }
    pGeantTrack->setDebugPrintFlag(debugPrintFlag);
    nGSegments += pGeantTrack->testGeantTrack(pGeantKine);
    geantBugFlag = pGeantTrack->isGeantBug();
    if(debugPrintFlag) pGeantTrack->print();
    if(!pGeantTrack->isInMdc()) continue;
    if(pGeantTrack->getNWires()==0) continue;
    nTracks++;
  }
  return kTRUE;
}

void HMdcGeantEvent::print(void) {
  for(Int_t i=0;i<nTracks;i++) ((HMdcGeantTrack*)At(i))->print();
}

void HMdcGeantEvent::clearOSegClus(void) {
  // This method is used for outer clusters when magnet on
  // for each inner segment.
  for(Int_t ind=0;ind<sizeGSegArr;ind++) {
    HMdcGeantSeg* pMdcGSeg = (HMdcGeantSeg*)(geantSegments.At(ind));
    if(pMdcGSeg==NULL || pMdcGSeg->getIOSeg()!=1) continue;
    pMdcGSeg->clearClus();
    pMdcGSeg->clearSeg();
  }
}

void HMdcGeantEvent::clearOClus(void) {
  // This method is used for outer clusters when magnet on
  // for each inner segment.
  for(Int_t ind=0;ind<sizeGSegArr;ind++) {
    HMdcGeantSeg* pMdcGSeg = (HMdcGeantSeg*)(geantSegments.At(ind));
    if(pMdcGSeg==NULL || pMdcGSeg->getIOSeg()!=1) continue;
    pMdcGSeg->clearClus();
  }
}

void HMdcGeantEvent::clearOSeg(void) {
  // This method is used for outer clusters when magnet on
  // for each inner segment.
  for(Int_t ind=0;ind<sizeGSegArr;ind++) {
    HMdcGeantSeg* pMdcGSeg = (HMdcGeantSeg*)(geantSegments.At(ind));
    if(pMdcGSeg==NULL || pMdcGSeg->getIOSeg()!=1) continue;
    pMdcGSeg->clearSeg();
  }
}

void HMdcGeantEvent::getEventVertex(Float_t &ax, Float_t &ay, Float_t &az) const {
  ax = xVertex;
  ay = yVertex;
  az = zVertex;
}

//===========================================================================

void HMdcTrackInfSim::fillClusTrackInf(HMdcClusSim* pClusSim) {
  sector        = pClusSim->getSec();
  segment       = pClusSim->getIOSeg();
  nWiresTot     = pClusSim->getNDrTimes();
  listTmp       = *((HMdcList12GroupCells*)pClusSim);
  modInd        = -2;
  isWrCollected = kFALSE;
  Int_t indPar  = pClusSim->getIndexParent();
  HMdcClusSim*  pParentClus = NULL;
  if(indPar>=0) pParentClus = (HMdcClusSim*)HMdcGetContainers::getObject()->
                              getCatMdcClus()->getObject(indPar);

  collectTracksInf();

  if( !isWrCollected ) pClusSim->calcTrList(); // old calculation
  else {
    Int_t        modICl  = pClusSim->getMod();
    if(modICl>0) modICl &= 1;
    Int_t typeClFinder = pClusSim->getTypeClFinder();
    Int_t maxNTrcks    = pClusSim->getArrsSize();
    Int_t nTrcks       = (numTracks < maxNTrcks) ? numTracks : maxNTrcks;
    for(Int_t n=0; n<nTrcks; n++) {
      Int_t ind = sortedInd[n];
      Int_t bin = pClusSim->addTrack(tracksNum[ind],numWires[ind],
          wiresList[ind].getListLayers(0),wiresList[ind].getListLayers(1));
      HMdcGeantSeg* gntSeg = gntSegList[ind];
      if(gntSeg != NULL) {
        if(typeClFinder == 1 && modICl >= 0)
             pClusSim->setNDigiTimes(bin,gntSeg->getNDrTimesMod(modICl));
        else pClusSim->setNDigiTimes(bin,gntSeg->getNDrTimes());

        gntSeg->analyseClust(pClusSim,bin,modICl);
        if(pParentClus) testIOMatching(ind,pClusSim,bin,pParentClus);
      } else {
        pClusSim->setNDigiTimes(bin,0);
        UChar_t status = ind == embedInd ? DEFAULT_SET : 0;
        if(pParentClus && ind == embedInd) {
          Int_t embRealTrack   = gHades->getEmbeddingRealTrackId();
          Int_t indPar         = pParentClus->getTrackIndex(embRealTrack);
          if(indPar >= 0) {
            status |= SEGS_MATH_OK;
            pParentClus->setTrackStatus(indPar, pParentClus->getTrackStatus(indPar) | SEGS_MATH_OK);
          }
        }
        pClusSim->setTrackStatus(bin,status);
      }
    }
  }
  pClusSim->cleanRest();
}

void HMdcTrackInfSim::testIOMatching(Int_t ind,HMdcClusSim* pClusSim,Int_t bin,
    HMdcClusSim* pClusSimPar) {
  // Check of inner-outer segments matching:
  // ind - segment index in arrays of HMdcTrackInfSim
  // bin - segment index in HMdcClusSim tracks list
  if(pClusSimPar == 0) return;
  HMdcGeantSeg*   gntOuterSeg = gntSegList[ind];
  HMdcGeantTrack* gntTrack    = gntTrackList[ind];
  Int_t  track      = gntTrack->getTrack();
  Int_t  nTrPar     = pClusSimPar->getNTracks();
  Bool_t isMatching = kFALSE;
  Int_t binPar = -1;
  if(!pClusSim->isFakeContribution(bin)) {
    for(Int_t trInd=0;trInd<nTrPar;trInd++) {
      if(pClusSimPar->getTrack(trInd) != track) continue;
//---???      if(pClusSimPar->isFakeContribution(trInd)) continue;

      Char_t dir   = gntOuterSeg->getDirection();
      Int_t segInd = segIngGTrack[ind];
      if(dir==1) segInd--; // Tr.dir.: innerSegment->outerSeg
      else       segInd++; // Tr.dir.: outerSeg->innerSegment
      HMdcGeantSeg* pInnerGSeg = gntTrack->getSegment(segInd);
      if(pInnerGSeg == 0) continue;
      if(pInnerGSeg->getIOSeg() != 0) continue;
      if(dir != pInnerGSeg->getDirection()) continue;
      if(pInnerGSeg->nIdentDrTimes(pClusSimPar)<2) continue;
      binPar = trInd;
      isMatching=kTRUE;
      break;
    }
  }
  if( isMatching ) {
    pClusSim->setTrackStatus(bin,pClusSim->getTrackStatus(bin) | SEGS_MATH_OK);
    pClusSimPar->setTrackStatus(binPar,pClusSimPar->getTrackStatus(binPar) | SEGS_MATH_OK);
  }  
}

void HMdcTrackInfSim::testIOMatching(Int_t ind,HMdcSegSim* pSegSim,Int_t bin,HMdcSegSim* pSegPar) {
  // Check of inner-outer segments matching:
  // ind - segment index in arrays of HMdcTrackInfSim
  // bin - segment index in pSegSim tracks list
  // pSegPar - parent segment for the pSegSim
  if(pSegPar == 0) return;
  HMdcGeantSeg*   gntOuterSeg = gntSegList[ind];
  HMdcGeantTrack* gntTrack    = gntTrackList[ind];
  Int_t  track      = gntTrack->getTrack();
  Int_t  nTrPar     = pSegPar->getNTracks();
  Bool_t isMatching = kFALSE;
  Int_t binPar = -1;
  HMdcCellGroup12 geantGr;
  if(!pSegSim->isFakeContribution(bin)) {
    for(Int_t trInd=0;trInd<nTrPar;trInd++) {
      if(pSegPar->getTrack(trInd) != track) continue;
//---???      if(pSegPar->isFakeContribution(trInd)) continue;
      Char_t dir   = gntOuterSeg->getDirection();
      Int_t segInd = segIngGTrack[ind];
      if(dir==1) segInd--; // Tr.dir.: innerSegment->outerSeg
      else       segInd++; // Tr.dir.: outerSeg->innerSegment
      HMdcGeantSeg* pInnerGSeg = gntTrack->getSegment(segInd);
      if(pInnerGSeg == 0) continue;
      if(pInnerGSeg->getIOSeg() != 0) continue;
      if(dir != pInnerGSeg->getDirection()) continue;
      for(Int_t l=0; l<12; l++) geantGr.setLayerGroup(l,pInnerGSeg->getOneLayerGroup(l));
      if(geantGr.getNSharedCells(pSegPar)<2) continue;
      
      binPar = trInd;
      isMatching=kTRUE;
      break;
    }
  }
  if( isMatching ) {
    pSegSim->setTrackStatus(bin,pSegSim->getTrackStatus(bin) | SEGS_MATH_OK);
    pSegPar->setTrackStatus(binPar,pSegPar->getTrackStatus(binPar) | SEGS_MATH_OK);
  }  
}

// ??????????????? Is it needed ???????????????????:
void HMdcTrackInfSim::addClusModTrackInf(HMdcClusSim* pClusSim) {
  // This function can be used ONLY after
  // fillClusTrackInf(HMdcClusSim* pClusSim) for the same pClusSim!!!
  if( !isWrCollected) {
    pClusSim->calcTrListMod(*pClusSim,0); // old calculation
    pClusSim->calcTrListMod(*pClusSim,1); // old calculation
  } else {
    Short_t nWires[numTracks];
    Int_t sInd[numTracks];
    for(Int_t m=0;m<2;m++) {
      for(Int_t tr=0;tr<numTracks;tr++) nWires[tr] = wiresList[tr].getNDrTimesMod(m);
      TMath::Sort(numTracks, nWires, sInd);
      Int_t maxTrcks     = pClusSim->getArrsSize();
      Int_t nSettedTrcks = 0;
      for(Int_t n=0; n<numTracks && nSettedTrcks<maxTrcks; n++) {
        Int_t ind = sInd[n];
        if(nWires[ind] == 0) break;
        pClusSim->setTrackM(m,nSettedTrcks,tracksNum[ind],nWires[ind]);
        nSettedTrcks++;
      }
      pClusSim->setNTracksM(m,nSettedTrcks);
    }
  }
}

void HMdcTrackInfSim::fillClusModTrackInf(HMdcClusSim* pClusSim,
    HMdcList12GroupCells* wrLst, Int_t modi) {
  // sector and segment in pClusSim must be setted already!
  isWrCollected = kFALSE;
  if(modi<0 || modi>1) return;
  sector    = pClusSim->getSec();
  segment   = pClusSim->getIOSeg();
  nWiresTot = pClusSim->getNDrTimesMod(modi);
  listTmp   = *wrLst;
  modInd    = modi;

  collectTracksInf();

  if( !isWrCollected ) pClusSim->calcTrListMod(*wrLst,modInd); // old calc.
  else {
    Int_t maxNTrcks    = pClusSim->getArrsSize();
    Int_t nSettedTrcks = (numTracks<maxNTrcks) ? numTracks : maxNTrcks;
    for(Int_t n=0; n<nSettedTrcks; n++) {
      Int_t ind = sortedInd[n];
      pClusSim->setTrackM(modInd,n,tracksNum[ind],numWires[ind]);
    }
    pClusSim->setNTracksM(modInd,nSettedTrcks);
  }
}

void HMdcTrackInfSim::fillSegTrackInf(HMdcSegSim* pSegSim,HMdcList24GroupCells* wrLst,
                                      HMdcSegSim* pSegPar) {
  isWrCollected = kFALSE;
  sector        = pSegSim->getSec();
  segment       = pSegSim->getIOSeg();
  nWiresTot     = pSegSim->getSumWires();
  if(!wrLst->getSeg(listTmp,segment)) return;
  modInd        = -2;

  collectTracksInf();

  if( !isWrCollected ) pSegSim->calcNTracks(); // old calculation
  else {
    if(pSegSim->getNTracks() != 0) pSegSim->clearSimInfo();
    Int_t   nSettedTrcks = numTracks <= 5 ? numTracks : 5;
    for(Int_t n=0; n<nSettedTrcks; n++) {
      Int_t   ind    = sortedInd[n];
      UChar_t nWires = numWires[ind] < 256 ? numWires[ind] : 255;
      HMdcGeantSeg* gntSeg = gntSegList[ind];
      if(gntSeg != NULL) {
        pSegSim->addTrack(tracksNum[ind],nWires,gntSeg->getNDrTimes());
        gntSeg->analyseSeg(pSegSim,n,wiresList[ind],
                           pWiresArr->calcChi2(wiresList[ind]),nWiresTot-nWires);
// No IOMatch        if(pSegPar) testIOMatching(ind,pSegSim,n,pSegPar);
      } else {
        pSegSim->addTrack(tracksNum[ind],nWires,0);
        UChar_t status = ind == embedInd ? DEFAULT_SET : SEGS_MATH_OK; //0; No match. info any more!
//         if(pSegPar && ind == embedInd) {
//           Int_t embRealTrack   = gHades->getEmbeddingRealTrackId();
//           Int_t indPar         = pSegPar->getTrackIndex(embRealTrack);
//           if(indPar >= 0) {
//             status |= SEGS_MATH_OK;
//             pSegPar->setTrackStatus(indPar, pSegPar->getTrackStatus(indPar) | SEGS_MATH_OK);
//           }
//         }
        pSegSim->setTrackStatus(n,status);
      }
    }
  }
}

void HMdcTrackInfSim::fillHitTrackInf(HMdcHitSim* pHitSim,HMdcList24GroupCells* wrLst) {
  isWrCollected = kFALSE;
  sector        = pHitSim->getSector();
  modInd        = pHitSim->getModule() & 1;
  segment       = pHitSim->getModule() >> 1;
  nWiresTot     = pHitSim->getSumWires();
  if(!wrLst->getSeg(listTmp,segment)) return;

  collectTracksInf();

  if( !isWrCollected ) pHitSim->calcNTracks(); // old calculation
  else {
    Int_t   segListTr[numTracks];
    UChar_t segNTimes[numTracks];
    for(Int_t n=0; n<numTracks; n++) {
      Int_t ind    = sortedInd[n];
      segListTr[n] = tracksNum[ind];
      segNTimes[n] = (numWires[ind] < 256) ? numWires[ind] : 255;
    }
    pHitSim->setNTracks(numTracks,segListTr,segNTimes);
  }
}

Bool_t HMdcTrackInfSim::collectTracksInf(UChar_t sec,UChar_t seg, HMdcList12GroupCells* wrLst) {
  sector        = sec;
  segment       = seg;
  listTmp       = *wrLst;
  modInd        = -2;
  isWrCollected = kFALSE;

  if(!collectTracksInf()) return kFALSE;
  return kTRUE;
}

Bool_t HMdcTrackInfSim::collectTracksInf(void) {
  HMdcEvntListCells* pMdcListCells = HMdcEvntListCells::getExObject();
  if(pMdcListCells == 0) return kFALSE;
  HMdcGeantEvent* pGeantEvent = HMdcGeantEvent::getExObject();
  if(pGeantEvent == 0) return kFALSE;

  HMdcSecListCells& fSecListCells = (*pMdcListCells)[sector];
  Int_t embRealTrack   = gHades->getEmbeddingRealTrackId();
  numTracks = 0;
  embedInd  = -1;
  noiseInd  = -1;
//----------  pWiresArr = NULL;

  Int_t layStr = (modInd != 1) ?  0 : 6;
  Int_t layEnd = (modInd != 0) ? 12 : 6;
  for(Int_t lay=layStr; lay<layEnd; lay++) {
    Int_t   mod  = lay/6 + segment*2;
    Int_t   cell = -1;
    UChar_t nTms = 0;
    HMdcLayListCells& fLayListCells = fSecListCells[mod][lay%6];
    while((cell=listTmp.next(lay,cell,nTms)) >= 0) {
      for(Int_t tm=1; tm<3; tm++) {
        Int_t track = fLayListCells.getGeantTrack(cell,tm);
        if(track == 0) continue;
        if(track<0 || track == embRealTrack) {
          // Embeded real track and noise --------------------------------------
          if(track != embRealTrack) track = -99;
          listTmp.delTime(lay,cell,tm);
          Int_t& ind = (track == embRealTrack) ? embedInd : noiseInd;
          if(ind < 0) {
            if(numTracks < aSize) ind = numTracks++;
            else {
              for(Int_t i=0;i<aSize;i++) if(i!=noiseInd && i!=embedInd) {
                if(ind<0 || numWires[i]<numWires[ind]) ind = i;
              }
            }
            wiresList[ind].clear();
            tracksNum[ind]    = track;
            gntSegList[ind]   = NULL;
            gntTrackList[ind] = NULL;
            numWires[ind]     = 1;
            segIngGTrack[ind] = -1;
          } else numWires[ind]++;
          wiresList[ind].setTime(lay,cell,tm);
        } else {
          // Geant track -------------------------------------------------------
          HMdcGeantTrack* pGeantTrack = pGeantEvent->getGeantTrack(track);
          if(pGeantTrack == NULL) continue; // Geant bug;
          Short_t nGeantSegs = pGeantTrack->getNSegments();
          for(Int_t segn=0;segn<nGeantSegs;segn++) {
            HMdcGeantSeg* pGeantSeg = pGeantTrack->getSegment(segn);
            if( pGeantSeg->getSec()   != sector || pGeantSeg->getIOSeg() != segment ) continue;
            numWires[numTracks] = listTmp.compareAndUnset(pGeantSeg,wiresList+numTracks,modInd);
            if(numWires[numTracks] == 0) continue;
            Int_t ind = numTracks;
            if(numTracks < aSize) numTracks++;
            else {
              for(Int_t i=0;i<aSize;i++) if(i!=noiseInd && i!=embedInd) {
                if(numWires[i] < numWires[ind]) ind = i;
              }
              if(ind == aSize) continue;
              wiresList[ind]  = wiresList[numTracks];
              numWires[ind]   = numWires[numTracks];
            }
            tracksNum[ind]    = track;
            gntSegList[ind]   = pGeantSeg;
            gntTrackList[ind] = pGeantTrack;
            segIngGTrack[ind] = segn;
          }
        }
      }
    }
  }

  isGntBugEvent       = pGeantEvent->isGeantBug();
  numWires[numTracks] = listTmp.getNDrTimes(layStr,layEnd-1);
  if(numWires[numTracks] == 0) isGntBugSeg = kFALSE;
  else {
    HLocation loc;
    loc.set(4,sector,0,0,0);
    HCategory* cat = HMdcGetContainers::getObject()->getCatMdcCal1();
    Int_t layInd=0;
    Int_t cell=-1;
    while(listTmp.getNextCell(layInd,cell)) {
      UChar_t times = listTmp.getTime(layInd,cell);
      loc.set(4,sector,layInd/6+segment*2,layInd%6,cell);
      HMdcCal1Sim* pMdcCal1Sim = (HMdcCal1Sim*)cat->getObject(loc);
      Int_t tr1 = (times&1)!=0 ? pMdcCal1Sim->getNTrack1() : 0;
      Int_t tr2 = (times&2)!=0 ? pMdcCal1Sim->getNTrack2() : 0;
      for(Int_t ind=0;ind<numTracks;ind++) {
        if(tracksNum[ind] == tr1) {
          numWires[ind]++;
          numWires[numTracks]--;
          listTmp.delTime(layInd,cell,1);
          tr1 = -1;
        }
        if(tracksNum[ind] == tr2) {
          numWires[ind]++;
          numWires[numTracks]--;
          listTmp.delTime(layInd,cell,2);
          tr2 = -1;
        }
      }
      if(numWires[numTracks] == 0) break;
//       if(tr1 > 0 && numTracks < aSize) {
//         Int_t ind = numTracks;
//         numTracks++; 
//         numWires[numTracks] = numWires[ind];
//         wiresList[ind].clear();
//         numWires[ind]  = 1;
//         tracksNum[ind] = tr1;
//         numWires[numTracks]--;
//         wiresList[ind].setTime(layInd,cell,1);
//         listTmp.delTime(layInd,cell,1);
//         gntSegList[ind]   = NULL;
//         gntTrackList[ind] = NULL;
//         segIngGTrack[ind] = -1;
//       }
//       if(tr2 > 0 && numTracks < aSize) {
//         Int_t ind = numTracks;
//         numTracks++; 
//         numWires[numTracks] = numWires[ind];
//         wiresList[ind].clear();
//         numWires[ind]  = 1;
//         tracksNum[ind] = tr2;
//         numWires[numTracks]--;
//         wiresList[ind].setTime(layInd,cell,2);
//         listTmp.delTime(layInd,cell,2);
//         gntSegList[ind]   = NULL;
//         gntTrackList[ind] = NULL;
//         segIngGTrack[ind] = -1;
//       }
    }
  }
    
isGntBugSeg = kFALSE;  // The case NumWires/layer > 48 ignored now partially !!! 
//   if(numWires[numTracks] == 0) isGntBugSeg = kFALSE;
//   else {
//     
//     isGntBugSeg = kTRUE;
//     Int_t ind = numTracks;
//     if(numTracks < aSize) numTracks++;
//     else {
//       for(Int_t i=0;i<aSize;i++) if(i!=noiseInd && numWires[i]<numWires[ind]) ind = i;
//       if(ind == embedInd) embedInd = -1;
//       numWires[ind]  = numWires[numTracks];
//     }
//     tracksNum[ind]    = -9;  // GEANT bug!
// printf("-9 => %i  nw=%i\n",tracksNum[ind],numWires[ind]);
// // Int_t lind=0;
// // Int_t c=-1;
// //  HCategory* cat = HMdcGetContainers::getObject()->getCatMdcCal1();
// //   HLocation loc;
// //   loc.set(4,sector,0,0,0);
// // while(listTmp.getNextCell(lind,c)) {
// //     loc[1]=lind/6+segment*2;
// //     loc[2]=lind%6;
// //     loc[3]=c;
// //  HMdcCal1Sim* cal1=(HMdcCal1Sim*)cat->getObject(loc);
// // printf("%i %i %i %i: tttt=%i tof=%f\n",sector,loc[1],lind%6,c,cal1->getNTrack1(),cal1->getTof1());
// // }
//     wiresList[ind]    = listTmp;
//     gntSegList[ind]   = NULL;
//     gntTrackList[ind] = NULL;
//   }

  TMath::Sort(numTracks, numWires, sortedInd);
  isWrCollected = kTRUE;
  return kTRUE;
}

Int_t HMdcTrackInfSim::getTrack(Int_t trInd) const {
  if(trInd>=0 && trInd<numTracks) return 0;
  return tracksNum[sortedInd[trInd]];
}

Short_t HMdcTrackInfSim::getNumWires(Int_t trInd) const {
  if(trInd>=0 && trInd<numTracks) return 0;
  return numWires[sortedInd[trInd]];
}

UChar_t HMdcTrackInfSim::getNumLayers(Int_t trInd, Int_t modi) const {
  if(trInd>=0 && trInd<numTracks) return 0;
  Int_t indx = sortedInd[trInd];
  if(modi==0 || modi==1) return HMdcBArray::getNSet(wiresList[indx].getListLayers(modi));
  return HMdcBArray::getNSet(wiresList[indx].getListLayers(0))+
         HMdcBArray::getNSet(wiresList[indx].getListLayers(1));
}

HMdcList12GroupCells* HMdcTrackInfSim::getListCells(Int_t trInd) {
  if(trInd>=0 && trInd<numTracks) return 0;
  return &(wiresList[sortedInd[trInd]]);
}

HMdcGeantSeg* HMdcTrackInfSim::getMdcGeantSeg(Int_t trInd) {
  if(trInd>=0 && trInd<numTracks) return 0;
  return gntSegList[sortedInd[trInd]];
}

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