//Task for all Matchings, fills hmetamatch
//*-- Author : Anar Rustamov (27.08.2003)
//*-- Modified : 14/10/2004 by V. Pechenov
//*-- Modified : 07/12/2007 by A. Rustamov
using namespace std;

#include "hmetamatchF.h"
#include "hevent.h"
#include "hcategory.h"
#include "hiterator.h"
#include "hspectrometer.h"
#include "hruntimedb.h"
#include "hmdcdetector.h"
#include "hspecgeompar.h"
#include "tofdef.h"
#include "richdef.h"
#include "rpcdef.h"
#include "hmatrixcategory.h"
#include "hbasetrack.h"
#include "hshowergeometry.h"
#include "hmetamatch.h"
#include "hgeomtransform.h"
#include "hmdctrackddef.h"
#include "hmdctrackgdef.h"
#include "hmdcseg.h"
#include "hmdctrkcand.h"
#include "hmdctrackgspline.h"
#include "hmdctrackgcorrpar.h"
#include "hgeomvector.h"
#include "hlocation.h"
#include "hshowerhittoftrack.h"
#include "hmdcgetcontainers.h"
#include "htofhit.h"
#include "htofcluster.h"
#include "htofgeompar.h"
#include "hrpcgeompar.h"
#include "hrichhit.h"
#include "hmetamatchpar.h"
#include "hmdcsizescells.h"
#include <iostream>

ClassImp(HMetaMatchF)

HMetaMatchF::HMetaMatchF() {
  setInitParam();
}

HMetaMatchF::HMetaMatchF(const Text_t *name,const Text_t *title):
  HReconstructor(name,title) {
  setInitParam();
}

void HMetaMatchF::setInitParam(void) {
  fMatchPar      = NULL;
  fCatTrkCand    = NULL;
  fCatMdcSeg     = NULL;
  fCatRich       = NULL;
  Spline         = NULL;
  fTrkCandIter   = NULL;
  fTofGeometry   = NULL;
  IterTof        = NULL;
  IterTofCluster = NULL;
  IterRich       = NULL;
  fTrkCandIter   = NULL;
  IterShower     = NULL;
  fShrGeometry   = NULL;
  fRpcGeometry   = NULL;
  fCatRpcCluster = NULL;
  IterRpcCluster = NULL;
  stNotMatTracks = kFALSE;
  for(Int_t s=0;s<6;s++) {labTrans[s]=0;}
  Spline=0;
  rad2deg = TMath::RadToDeg();
  locTofHit.set(3,0,0,0);
}

HMetaMatchF::~HMetaMatchF() {
     HMdcSizesCells::deleteCont();
   //HMdcGetContainers::deleteCont();

   if(fTrkCandIter) {
      delete fTrkCandIter;
      fTrkCandIter=0;
   }
  if(IterRich){
    delete IterRich;
    IterRich=0;
  }
  if(IterTof) {
    delete IterTof;
    IterTof=0;
  }
  if(IterTofCluster) {
    delete IterTofCluster;
    IterTofCluster=0;
  }
  if(IterShower) {
    delete IterShower;
    IterShower=0;
  }
}

Bool_t HMetaMatchF::init() {
  if (!gHades) return kFALSE;
  HRuntimeDb *rtdb=gHades->getRuntimeDb();
  if(!rtdb) return kFALSE;
  HEvent *event=gHades->getCurrentEvent();
  if(!event) return kFALSE;
  
  HMdcTrackGCorrPar* corr=(HMdcTrackGCorrPar*)(rtdb->getContainer("MdcTrackGCorrPar"));
  if(corr){
      Spline=corr->getSPline();
  } else {
      Error("init()","ZERO pointer for MdcTrackGCorrPar received!");
      return kFALSE;
  }
  fGetCont=HMdcGetContainers::getObject();
  pSizesCells = HMdcSizesCells::getObject();
  fGetCont->getMdcGeomPar();  // for spline ???
  fGetCont->getSpecGeomPar(); // for lab-sec. transf
  //Spline=new HMdcTrackGSpline("Spline","Spline");
  fCatTrkCand=event->getCategory(catMdcTrkCand);
  if (!fCatTrkCand) {
     Error("init","NO catMdcTrkCand in input! STOP!!!");
     return kFALSE;
  }
  fTrkCandIter=(HIterator*)fCatTrkCand->MakeIterator();
  if(!fTrkCandIter) return kFALSE;
  fCatMdcSeg=event->getCategory(catMdcSeg);

  // TOF:
  fCatTof=event->getCategory(catTofHit);
  if(!fCatTof) Warning("init",
      "No catTofHit in input! \n Matching with TofHits will be skipped!");
  else IterTof=(HIterator*)fCatTof->MakeIterator();
  fCatTofCluster=event->getCategory(catTofCluster);
  if(!fCatTofCluster) {
    Warning("init","NO catTofCluster in input! \n Matching with TofClusters will be skipped!");
  } else IterTofCluster=(HIterator*)fCatTofCluster->MakeIterator();
  if(fCatTof || fCatTofCluster) fTofGeometry=
      (HTofGeomPar *)rtdb->getContainer("TofGeomPar");

  // Shower:
  HCategory *fCatKine=event->getCategory(catGeantKine);
  if(!fCatKine) {
    fCatShower = event->getCategory(catShowerHitTof);
    if(!fCatShower) Warning("init",
        "NO catShowerHitTof in input! \n Matching with Shower will be skipped!");
  } else {
    fCatShower = event->getCategory(catShowerHitTofTrack);
    if(!fCatShower) Warning("init",
        "NO catShowerHitTofTrack in input! \n Matching with Shower will be skipped!");
  }
  if(fCatShower) {
    IterShower=(HIterator*)fCatShower->MakeIterator();
    fShrGeometry=(HShowerGeometry*)rtdb->getContainer("ShowerGeometry");
  }
 
  // RICH:
  fCatRich=event->getCategory(catRichHit);
  if(fCatRich) IterRich=(HIterator*)fCatRich->MakeIterator();
  else Warning("HMetaMatchF::init()","NO RICH catRichHit in input! \n Matching with Rich will be skipped!");

  fCatMetaMatch=event->getCategory(catMetaMatch);
  if(!fCatMetaMatch) {
    Int_t size[2]={6,8000};
    fCatMetaMatch=new HMatrixCategory("HMetaMatch",2,size,0.5);
    if(fCatMetaMatch)
       event->addCategory(catMetaMatch,fCatMetaMatch,"Tracks");
  }
  fMatchPar=(HMetaMatchPar*)rtdb->getContainer("MetaMatchPar");

  return kTRUE;
}

Bool_t HMetaMatchF::reinit() {

  //Spline=corr->getSPline();
  if(!pSizesCells->initContainer()) return kFALSE;
  if (pSizesCells->hasChanged()) {
     // Geometry transformation from module to sector coord system for MDCs
     for(Int_t is=0; is<6; is++) {
	HMdcSizesCellsSec& fSCSec = (*pSizesCells)[is];
	if(&fSCSec == 0) continue; // sector is not active
	//secTrans[is]=*(fSCSec.getLabTrans());
	for(Int_t im=0; im<4; im++) {
	   HMdcSizesCellsMod& fSCMod=fSCSec[im];
	   if (&fSCMod) {
	       const HGeomTransform* tr = fSCMod.getSecTrans();
	       Spline->takeMiddleParams(tr,is,im);
	   }
	}
     }
  }
  Spline->initMiddleParamsAll();
  Spline->setKickPointer(&kickplane);
  for(Int_t s=0; s<6; s++) {
     if(fGetCont->getMdcDetector()->isSectorActive(s)) {
	labTrans[s]=&(fGetCont->getLabTransSec(s));
	if(labTrans[s]==0) return kFALSE;
     } else labTrans[s]=0;
  }
  if(fTofGeometry && !HMdcGetContainers::isInited(fTofGeometry)) return kFALSE;
  if(fShrGeometry && !HMdcGetContainers::isInited(fShrGeometry)) return kFALSE;

  if(fMatchPar==0 || !HMdcGetContainers::isInited(fMatchPar)) {
    Error("reinit","no parameters for matching!");
    return kFALSE;
  }
  setMatchingParam();
  return kTRUE;
}

void HMetaMatchF::setMatchingParam(void) {
  // Matching parameters.
  for(Int_t s=0;s<6;s++) {
    dThRich[s]=fabs( fMatchPar->getRichThetaMaxCut(s) -
                     fMatchPar->getRichThetaMinCut(s) )/2.;
    dPhRich[s]=fMatchPar->getRichSigmaPhi(s);
    dPhRichOff[s]=fMatchPar->getRichSigmaPhiOffset(s);
    qualityRichCut[s]=fMatchPar->getRichQualityCut(s);

    sigma2MdcInShrX[s]=
        fMatchPar->getShowerSigmaXMdc(s)*fMatchPar->getShowerSigmaXMdc(s);
    sigma2MdcInShrY[s]=
        fMatchPar->getShowerSigmaYMdc(s)*fMatchPar->getShowerSigmaYMdc(s);
    sShowerX[s]=fMatchPar->getShowerSigmaXOffset(s);
    sShowerY[s]=fMatchPar->getShowerSigmaYOffset(s);
    qualitySHOWERCut[s]=fMatchPar->getShowerQualityCut(s);

    sigmaTofX[s]=fMatchPar->getTofSigmaX(s);
    sigmaTofY[s]=fMatchPar->getTofSigmaY(s);
    sTofX[s]=fMatchPar->getTofSigmaXOffset(s);
    sTofY[s]=fMatchPar->getTofSigmaYOffset(s);
    qualityTOFCut[s]=fMatchPar->getTofQualityCut(s);
    richThetaMinCut[s]=fMatchPar->getRichThetaMinCut(s);
    richThetaMaxCut[s]=fMatchPar->getRichThetaMaxCut(s);
  }
}

Bool_t HMetaMatchF::finalize() {
  return kTRUE;
}

Int_t HMetaMatchF::execute() {
  if(fCatTrkCand->getEntries()<=0) return 0;

  collectShowerHits();
  collectTofHits();

  fTrkCandIter->Reset();
  HMdcTrkCand *pTrkCand=0;

  for(Int_t sec=0;sec<6;sec++) {
    setCurrentSector(sec);
    if(!fTrkCandIter->gotoLocation(sectorLoc)) continue;

    for(Int_t nTof=0; nTof<nTofHitsSec;nTof++) nTrCndForTHit[nTof]=0;
    for(Int_t nSh=0; nSh<nShowerHitsSec;nSh++) nTrCndForSHit[nSh]=0;
    Bool_t isCandInSec=kFALSE;
    while((pTrkCand=(HMdcTrkCand*)(fTrkCandIter->Next()))!=0) {
      if(pTrkCand->getNCandForSeg1()<0) continue;
      Short_t index1=pTrkCand->getSeg1Ind();
      if(index1<0) continue;
      segments[0]=(HMdcSeg*)fCatMdcSeg->getObject(index1);
      if(!segments[0]) continue;

      makeRichMatching();                 // MDC-RICH matching
      isCandInSec=kTRUE;
      if(!makeOuterSegMatch(pTrkCand))   //Outer segment matching 
      {
        break;
      }
      
    }
    if(isCandInSec) 
    {
      fillTofClstInd();
      fillMetaHitsInd();
    }
  }
  return 0;
}

Bool_t HMetaMatchF::makeOuterSegMatch(HMdcTrkCand* pTrkCand) {
  // Matching outer mdc segments with meta hits.
  // Inner mdc segment is the same.
  Bool_t isMMatch = kFALSE;
  for(Int_t nTof=0; nTof<nTofHitsSec;nTof++) {
    Int_t& nTrCndForT = nTrCndForTHit[nTof];
    if(nTrCndForT>=300) continue;
    mmTof[nTof][nTrCndForT]=-1;
    nTrCndForT++;
  }
  for(Int_t nSh=0; nSh<nShowerHitsSec;nSh++)  {
    Int_t& nTrCndForS=nTrCndForSHit[nSh];
    if(nTrCndForS>=300) continue;
    mmShr[nSh][nTrCndForS]=-1;
    nTrCndForS++;
  }
  

isPrint=kTRUE;
isPrint0=kTRUE;
  while((trkCandIndex=pTrkCand->getNextCandInd())>=0) {
    pTrkCand=(HMdcTrkCand*)fCatTrkCand->getObject(trkCandIndex);
    // MDC-Meta matching with outer mdc segment:
    meta = 0;
    Short_t index2=pTrkCand->getSeg2Ind();
    if(index2<0) continue;               // It can happen by error only.
    segments[1]=(HMdcSeg*)fCatMdcSeg->getObject(index2);
    if(!segments[1]) continue;           // It can happen by error only.
    if(segments[1]->getHitInd(1)<0)
    {
       Spline->calcSegPoints123(segments,mdcTrackPar);
       }
       else if(segments[1]->getHitInd(0)<0)
       {
       Spline->calcSegPoints123P4(segments,mdcTrackPar);
       }
    else
    {
       Spline->calcSegPoints(segments,mdcTrackPar);
    }
    for(Int_t p=0;p<4;p++) mdcTrackPar[p] *= 10.; // cm => mm
    mdcTrackPar[2] = secLabTrans->transFrom(mdcTrackPar[2]);
    mdcTrackPar[3] = secLabTrans->transFrom(mdcTrackPar[3]);

    Bool_t isNotMatch = kTRUE;
    for(Int_t nTof=0; nTof<nTofHitsSec;nTof++) {
      Int_t& nTrCndForT=nTrCndForTHit[nTof];
      if(nTrCndForT>=300) continue;

      //Variables to store x/y-deviation of META hit and segment normalized to META resolution
      qualTof=qualityTof(nTof,meta_dx,meta_dy);
      if(qualTof>qualityTOFCut[sector]) 
      {
      //if(mmTof[nTof][nTrCndForT-1]<0) mmTof[nTof][nTrCndForT-1]=-20; 
      continue;
      }
      if(mmTof[nTof][nTrCndForT-1] < 0) nTrCndForT--;
      testIndex=fillMeta(-1,nTof,meta);
      if(testIndex==-1) {isPrint=kFALSE;}
      mmTof[nTof][nTrCndForT]=testIndex;
      nTrCndForT++;
      isNotMatch = kFALSE;
      isMMatch   = kTRUE;
    }
    for(Int_t nSh=0; nSh<nShowerHitsSec;nSh++) {
      Int_t& nTrCndForS=nTrCndForSHit[nSh];
      if(nTrCndForS>=300) continue;
      qualShower=qualityShower(nSh,meta_dx,meta_dy);
      if(qualShower>qualitySHOWERCut[sector]) 
      {
     // if(mmTof[nSh][nTrCndForS-1]<0 ) mmTof[nSh][nTrCndForS-1]=-20; 
      continue;
      }
      if(mmShr[nSh][nTrCndForS-1] < 0) nTrCndForS--;
      testIndex=fillMeta(nSh,-1,meta);
      if(testIndex==-1) {isPrint=kFALSE;}
      mmShr[nSh][nTrCndForS]=testIndex;
      nTrCndForS++;
      isNotMatch = kFALSE;
      isMMatch   = kTRUE;
    }
    // Creating of MetaMatch object for MdcTrkCand which has not have matching to Meta
    if(isNotMatch && (nRichId>0 || stNotMatTracks)) {
      if(getMetaMatchSlot(0)==-1) isPrint=kFALSE;
      else isMMatch = kTRUE;
    }
  }
  
  // No outer segment and no meta hits:
  if(!isMMatch && meta==0) {
    if((getMetaMatchSlot(0))==-1) isPrint=kFALSE;
  }
  
  return kTRUE;
}

Int_t HMetaMatchF::getMetaMatchSlot(HMetaMatch* prevMM) {
  // For TrackCand with inner MDC segments only
  Int_t metaIndex=-1;
  HMetaMatch * metaNew=(HMetaMatch*)fCatMetaMatch->getNewSlot(sectorLoc,
      &metaIndex);
  if(!metaNew) { 
    if(isPrint){
    Warning("getMetaMatchSlot",
        "No slot available in sector %i. size of catMetaMatch is %i!",
        sectorLoc[0]+1,fCatMetaMatch->getEntries());
	}
    return -1;
  }
  if(prevMM==0) {
    HMdcTrkCand* pTrkCand=(HMdcTrkCand*)fCatTrkCand->getObject(trkCandIndex);
    if(pTrkCand) pTrkCand->setMetaMatchInd((Short_t)metaIndex);
    meta = new(metaNew) HMetaMatch(sector,trkCandIndex,metaIndex);
    if( nRichId > 0) {
      if(nRichId>3) nRichId=3;
      meta->setNCandForRich(nRichId);
      for(Int_t i = 0; i < nRichId; i++) meta->setARichInd(i,aRichIndTable[i]);
    }
  } else meta=new(metaNew) HMetaMatch(prevMM,metaIndex); // (meta,metaIndex) == meta->nextCandForMeta=metaIndex

  HMetaMatch *metaFirst=
      (HMetaMatch*)(fCatMetaMatch->getObject(meta->getFirstMMForSameTrCnd()));
  metaFirst->incrNumMMForSameTrkCnd();
  return metaIndex;
}

Double_t HMetaMatchF::qualityTof(Int_t hit,Float_t& dX_n, Float_t& dY_n) {
  HGeomTransform &tofModSys =
      fTofGeometry->getModule(sector,tofModuleSec[hit])->getLabTransform();
  HGeomVector p1=tofModSys.transTo(mdcTrackPar[2]);
  HGeomVector p2=tofModSys.transTo(mdcTrackPar[3]);
  xSegCrTof=(p1(2)*p2(0)-p1(0)*p2(2))/(p1(2)-p2(2));
  ySegCrTof=(p1(2)*p2(1)-p1(1)*p2(2))/(p1(2)-p2(2));
  Double_t Xtof=tofHitsSec[hit].getX();
  Double_t Ytof=tofHitsSec[hit].getY();
  Double_t dX=(Xtof-xSegCrTof-sTofX[sector])/sigmaTofX[sector];
  Double_t dY=(Ytof-ySegCrTof-sTofY[sector])/sigmaTofY[sector];
  dX_n=Xtof-xSegCrTof-sTofX[sector];
  dY_n=Ytof-ySegCrTof-sTofY[sector];
  return sqrt(dX*dX + dY*dY);
}

void HMetaMatchF::collectTofHits(void) {
  if(fCatTof) IterTof->Reset();
  HTofHit *pTofHit;
  for(Int_t s=0;s<6;s++) nTofHits[s]=0;
  if(!fCatTofCluster) {
    // No TofCluster category:
    if(fCatTof)
      while((pTofHit=(HTofHit*)(IterTof->Next()))!=0 ) addTofHit(pTofHit,0);
  } else {
    // TofCluster category exist:
    IterTofCluster->Reset();
    HTofCluster *pTofCluster;
    while((pTofCluster=(HTofCluster*)(IterTofCluster->Next()))!=0 ) {
      Int_t tofClSize=pTofCluster->getClusterSize();
      if(tofClSize>2) continue;                  // tofClSize<=2 only !
      addTofCluster(pTofCluster);
      if(tofClSize<2) continue;
      if(fCatTof==0) continue;
      Int_t sec = pTofCluster->getSector();
      Int_t mod=pTofCluster->getModule();
      Int_t cell=pTofCluster->getCell();
      while((pTofHit=(HTofHit*)(IterTof->Next()))!=0 ) {
        if(sec!=pTofHit->getSector() || mod!=pTofHit->getModule() ||
            cell!=pTofHit->getCell()) continue;
        if(tofClSize==2) {  // two TofHits adding, tofClSize==2 only !
          addTofHit(pTofHit,0);                  // first TofHit for TofCluster
          pTofHit=(HTofHit*)(IterTof->Next());   // second ...
          addTofHit(pTofHit,0);                  // clust.size==2 selected !!!
        } else {            // tofClSize>2 only HTofHits can be used
          addTofHit(pTofHit,-1);
          for(Int_t h=0;h<tofClSize-1;h++) {
            pTofHit=(HTofHit*)(IterTof->Next());
            addTofHit(pTofHit,-1);
          }
        }
        break;
      }
    }
  }
}

void HMetaMatchF::addTofCluster(HTofCluster* pTofCluster) {
  addTof(pTofCluster,fCatTofCluster->getIndex(pTofCluster),
      pTofCluster->getClusterSize());
}

void HMetaMatchF::addTofHit(HTofHit* pTofHit,Int_t clSize) {
  if(pTofHit==0) Error("addTofHit"," Pointer to HTofHit == 0 !");
  else addTof(pTofHit,fCatTof->getIndex(pTofHit),clSize);
}

void HMetaMatchF::addTof(HTofHit* pTofHit,Int_t index, Int_t clSize) {
  Int_t sec = pTofHit->getSector();
  Int_t& nTofHSec=nTofHits[sec];
  Float_t Xtof,Ytof,Ztof;
  pTofHit->getXYZLab(Xtof,Ytof,Ztof);
  HGeomVector& point=tofHits[sec][nTofHSec];
  point.setXYZ(Xtof,Ytof,Ztof);
  HModGeomPar *module=fTofGeometry->getModule(sec,pTofHit->getModule());
  if(module==0) {
    Error("addTof"," Can't get transformation for tof. %i sec. %imod",
        sec,pTofHit->getModule());
    return;
  }
  HGeomTransform &trans = module->getLabTransform();
  point=trans.transTo(point);    // transformation to coor.sys. of tof module

  indexTofHit[sec][nTofHSec]  = index;
  tofClustSize[sec][nTofHSec] = clSize;
  tofModule[sec][nTofHSec]    = pTofHit->getModule();
  nTofHSec++;
}

Double_t HMetaMatchF::qualityShower(Int_t hit,Float_t& dX_n, Float_t& dY_n)  {
  const HGeomTransform& showerModSys=
    fShrGeometry->getTransform(sector,fShowerHitsSec[hit]->getModule());
  HGeomVector p1 = showerModSys.transTo(mdcTrackPar[2]);
  HGeomVector p2 = showerModSys.transTo(mdcTrackPar[3]);
  
  Float_t Xshr, Yshr, Zshr; 
  fShowerHitsSec[hit]->getXY(&Xshr, &Yshr);
  /////////////////
  Zshr    = fShowerHitsSec[hit]->getZ();
  p1.Z() -= Zshr;
  p2.Z() -= Zshr;
  ////////////////
  xSegCrShr = (p1(2)*p2(0)-p1(0)*p2(2))/(p1(2)-p2(2));
  ySegCrShr = (p1(2)*p2(1)-p1(1)*p2(2))/(p1(2)-p2(2));
  Double_t dX       = Xshr-xSegCrShr-sShowerX[sector];
  Double_t dY       = Yshr-ySegCrShr-sShowerY[sector];
  Double_t dXsigma2 = fShowerHitsSec[hit]->getSigmaX();
  dXsigma2= dXsigma2*dXsigma2 + sigma2MdcInShrX[sector];
  Double_t dYsigma2 = fShowerHitsSec[hit]->getSigmaY();
  dYsigma2= dYsigma2*dYsigma2 + sigma2MdcInShrY[sector];
  dX_n = dX;
  dY_n = dY;
  return sqrt(dX*dX/dXsigma2 +dY*dY/dYsigma2);
}

void HMetaMatchF::collectShowerHits(void) {
  for(Int_t s=0;s<6;s++) nShowerHits[s]=0;
  if(!fCatShower) return;
  IterShower->Reset();
  HShowerHit *pShowerHit;
  while((pShowerHit = (HShowerHit*)(IterShower->Next()))!=0) {
    Int_t sec = pShowerHit->getSector();
    if(nShowerHits[sec]==200) continue; //!!!
    fShowerHits[sec][nShowerHits[sec]] = pShowerHit;
    indexShrHit[sec][nShowerHits[sec]] = fCatShower->getIndex(pShowerHit);
    nShowerHits[sec]++;
  }
}

void HMetaMatchF::makeRichMatching(void) {
  
  Float_t mdcPhi    = segments[0]->getPhi()*rad2deg + (sector!=5 ? sector*60.:-60.);
  Float_t mdcTheta  = segments[0]->getTheta()*rad2deg;
  Float_t sinMTheta = TMath::Sin(segments[0]->getTheta());
  
  // Matching with HRichHit's:
  nRichId = 0;
  if(fCatRich) {
    HRichHit* pRichHit;
    IterRich->Reset();
    while((pRichHit=(HRichHit*)(IterRich->Next()))!=0) {
      if(pRichHit->getSector() != sector) continue;
      Float_t dTheta  = pRichHit->getTheta()-mdcTheta;
      if(dTheta<richThetaMinCut[sector] || dTheta>richThetaMaxCut[sector]) continue;
      Float_t qualPhi = pRichHit->getPhi() - mdcPhi - dPhRichOff[sector];
      qualPhi = TMath::Abs(qualPhi *sinMTheta/dPhRich[sector]);
      if(qualPhi>qualityRichCut[sector]) continue;
      Short_t ind=fCatRich->getIndex(pRichHit);
      addRing(qualPhi,ind,aRichIndTable, nRichId);
    }
  }
}

void HMetaMatchF::addRing(Float_t quality, Short_t ind,
    Short_t* indTable,Int_t& nRich) {
  // Adding matched ring in sorted list
  if(nRich==0 || quality>=qualRich[nRich-1]) {
    // if ring is first in the list or ring must be added to the end of list:
    if(nRich>=RICH_TAB_SIZE) return;
    indTable[nRich] = ind;
    qualRich[nRich] = quality;
    nRich++;
  } else {
    // if ring must be inserted in list:
    for(Int_t i=0;i<nRich;i++) {
      if(quality>=qualRich[i]) continue;
      if(nRich<RICH_TAB_SIZE) nRich++;
      for(Int_t ish=nRich-1;ish>i;ish--) {  // Shift of rings in list
        indTable[ish]=indTable[ish-1];      // ...
        qualRich[ish]=qualRich[ish-1];      // ...
      }
      indTable[i]=ind;                      // Inserting ring in list
      qualRich[i]=quality;                  // ...
      return;
    }
  }
}

void HMetaMatchF::setCurrentSector(Int_t sec) {
  sector          = sec;
  sectorLoc.set(1,sector);
  secLabTrans     = labTrans[sector];

  nShowerHitsSec  = nShowerHits[sec];
  nTofHitsSec     = nTofHits[sec];

  indexShrHitSec  = indexShrHit[sec];
  fShowerHitsSec  = fShowerHits[sec];

  tofHitsSec      = tofHits[sec];
  indexTofHitSec  = indexTofHit[sec];
  tofModuleSec    = tofModule[sec];

  tofClustSizeSec = tofClustSize[sec];
}

Int_t HMetaMatchF::fillMeta(Int_t hShower, Int_t hTof, HMetaMatch* metaOld) {
  Int_t metaIndex = getMetaMatchSlot(metaOld);
  if(metaIndex>=0) {
    meta->setdxMeta(meta_dx);
    meta->setdyMeta(meta_dy);
    if(hShower>=0) {
      meta->setShowerHitInd(indexShrHitSec[hShower]);
      meta->setQualitySHOWER(qualShower);
      meta->setMdcSegCross(xSegCrShr,ySegCrShr);
    }
    if(hTof>=0) {
      meta->setTofHitInd(indexTofHitSec[hTof]);
      meta->setTofClusterSize(tofClustSizeSec[hTof]);
      meta->setQualityTOF(qualTof);
      meta->setMdcSegCross(xSegCrTof,ySegCrTof);
    }
  }
  return metaIndex;
}

void HMetaMatchF::fillTofClstInd(void) {
  for(Int_t nTof=0; nTof<nTofHitsSec;nTof++) {
    if(nTrCndForTHit[nTof]<1) continue;  // No mdc tracks for this tof cluster.
    if(tofClustSizeSec[nTof]<2) continue;
    Int_t mmClstInd = mmTof[nTof][0];
    if(mmClstInd<0) continue;
    HMetaMatch *pMetaMatch = 0;
    for(Int_t i=0;i<=tofClustSizeSec[nTof];i++) {
      if(nTof+i>=nTofHitsSec) {
        Error("fillTofClstInd","tof hit index out of range!");
        break;
      }
      Int_t mmInd=mmTof[nTof+i][0];
      if(mmInd<0) {
        if(isPrint0)
        Error("fillTofClstInd","mmInd <0 !");
        break;
      }
      if(pMetaMatch) pMetaMatch->setNextMMForTofHit(mmInd);
      pMetaMatch = getMMObj(mmInd);
      if(pMetaMatch==0) break;
      pMetaMatch->setIndMMForTofClst(mmClstInd);
    }
  }
}

void HMetaMatchF::fillMetaHitsInd(void) {
  for(Int_t nTof=0; nTof<nTofHitsSec;nTof++) {
    Int_t& nTrCndForT = nTrCndForTHit[nTof];
    if(nTrCndForT<1) continue;  // No mdc tracks for this tof hit.
    Int_t mmFirstInd = mmTof[nTof][0];
    if(mmFirstInd < 0) continue;
    HMetaMatch* pMMatch = 0;
    for(Int_t nTrCnd=0;nTrCnd<nTrCndForT;nTrCnd++) {
      pMMatch = setMetaHitsInd(pMMatch,mmFirstInd,mmTof[nTof][nTrCnd],nTrCndForT);
      if(pMMatch==0) break;
    }
  }
  for(Int_t nSh=0; nSh<nShowerHitsSec;nSh++) {
    Int_t& nTrCndForS = nTrCndForSHit[nSh];
    if(nTrCndForS<1) continue;  // No mdc tracks for this shower hit.
    Int_t mmFirstInd = mmShr[nSh][0];
    if(mmFirstInd < 0) continue;
    HMetaMatch* pMMatch=0;
    for(Int_t nTrCnd=0;nTrCnd<nTrCndForS;nTrCnd++) {
      pMMatch = setMetaHitsInd(pMMatch,mmFirstInd,mmShr[nSh][nTrCnd],nTrCndForS);
      if(pMMatch==0) break;
    }
  }
}

HMetaMatch* HMetaMatchF::setMetaHitsInd(HMetaMatch* pMMatch, Int_t mmFirstInd,
    Int_t mmInd, Int_t totNTrCnd) {
  // mmFirstInd<0 is tested before call setMetaHitsInd
  if(pMMatch) pMMatch->setNextMMForSameMeta(mmInd);
  if(mmInd<0) return 0;
  HMetaMatch* metaNext = getMMObj(mmInd);
  if(metaNext==0) return 0;
  if(pMMatch==0) metaNext->setNumMMForSameMeta(totNTrCnd);
  else        metaNext->setNumMMForSameMeta(-1);
  metaNext->setFirstMMForSameMeta(mmFirstInd);
  return metaNext;
}

HMetaMatch* HMetaMatchF::getMMObj(Int_t ind) {
  HMetaMatch* pMMatch=0;
  if(ind>=0) pMMatch=(HMetaMatch*)(fCatMetaMatch->getObject(ind));
  if(pMMatch==0) Error("getMMObj","no HMetaMatch object with index =%i !",ind);
  return pMMatch;
}