hparticlebtringf.cc

Go to the documentation of this file.

//*-- Author   : Patrick Sellheim, Georgy Kornakov
//*-- Created  : 09/01/2014

//_HADES_CLASS_DESCRIPTION
//////////////////////////////////////////////////////////////////////////////
//
//  HParticleBtRingF
//
//  This class predicts fired pads based on theta and phi angle of a track.
//  Predicted pads and rich cal hits are stored and used to search for
//  fired pads matching to prediction area. 
//
//
//
//
//////////////////////////////////////////////////////////////////////////////

#include "hparticlebtringf.h"
#include "TF2.h"
#include "hcategory.h"
#include "hparticlecand.h"
#include "heventheader.h"

#include "hrichcal.h"
#include "hades.h"
#include "hruntimedb.h"
#include "hcategorymanager.h"
#include "hparticlebtangletrafo.h"



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

ClassImp(HParticleBtRingF)

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

HParticleBtRingF::HParticleBtRingF( ){

}

HParticleBtRingF::~HParticleBtRingF( ){

}


Bool_t HParticleBtRingF::init() {

    fNSector = 6;

    fBtPar =(HParticleBtPar*) gHades->getRuntimeDb()->getContainer("ParticleBtPar");
    if(!fBtPar) {
        Error ("init()", "Retrieve ZERO pointer for ParticleBtPar!");
        return kFALSE;
    }
    //Init constants
    fNVertex                = fBtPar->getNVertex();
    fNParMean               = fBtPar->getNParMean();
    fNParSigma              = fBtPar->getNParSigma();
    fNRingSeg               = fBtPar->getNRingSegments();
    fNRingSegStep           = fBtPar->getRingSegStep();
    fNRingSegOffset         = fBtPar->getRingSegOffset();

    fNRichSeg               = fBtPar->getNRichSeg();
    fRichSegBorderX         = fBtPar->getRichSegBorderX();
    fRichSegBorderY         = fBtPar->getRichSegBorderY();
  
    fThetaAngleMin          = fBtPar->getThetaAngleMin();
    fThetaAngleMax          = fBtPar->getThetaAngleMax();
    fSizeMatrix             = fBtPar->getSizeMatrix();
    memcpy(&fChargeLimit[0]        ,fBtPar->getChargeThres()   , sizeof(fChargeLimit) );
    memcpy(&fChargeLimitMaximum[0] ,fBtPar->getChargeThresMax(), sizeof(fChargeLimitMaximum) );

    fNSigma                 = fBtPar->getSigmaRange();
    fMaxSigmaRange          = fBtPar->getSigmaRangeMax();

    fVertexPosMin           = fBtPar->getVertexPosMin();
    fVertexPosMax           = fBtPar->getVertexPosMax();
    fVertexStep             = fBtPar->getVertexStep();

    memcpy(&fPhiOff[0] , fBtPar->getPhiOffset() , sizeof(fPhiOff ) );
    memcpy(&fPhiOff2[0], fBtPar->getPhiOffset2(), sizeof(fPhiOff2 ) );
    fPhiOffsetPar        = fBtPar->getOffsetPar();


    fPol2DMean.resize(fNRingSeg);
    fPol2DSigma.resize(fNRingSeg);
    for(Int_t bin=0; bin<fNRingSeg; bin++){
        for(Int_t vertex=0; vertex<fNVertex; vertex++){
            //Theta,Phi Range, + x: x is half step width of binning used for calculating theses values
            fPol2DMean[bin].push_back( new TF2(Form("meanPol_%i_v%i",bin,vertex),"x*x*x*[0]+x*x*[1]+x*[2]+y*y*[3]+y*[4]+[5]+y*(x*x*[6]+x*[7])",20.+3.3,86.,0+7.51,60.)); 
            fPol2DSigma[bin].push_back(new TF2(Form("sigmaPol_%i_v%i",bin,vertex),"TMath::Power(x,3)*[0]+TMath::Power(x,2)*[1]+x*[2]+TMath::Power(y,4)*[3]+TMath::Power(y,3)*[4]+TMath::Power(y,2)*[5]+y*[6]+[7]+y*(x*x*[8]+x*[9]+[10]) + x*(y*y*y*[11]+y*y*[12])",20.+3.3,86,8+7.51,60.));
            for(Int_t par=0; par<fNParMean; par++)
                fPol2DMean[bin][vertex] ->SetParameter(par, fBtPar->getTF2ParMean( bin, vertex, par));
            for(Int_t par=0; par<fNParSigma; par++)
                fPol2DSigma[bin][vertex]->SetParameter(par, fBtPar->getTF2ParSigma(bin, vertex, par));
        }
        fRad2Deg.push_back(TMath::DegToRad()*(bin*fNRingSegStep+fNRingSegOffset));
        fRad2DegX.push_back(TMath::Sin(fRad2Deg[bin]));
        fRad2DegY.push_back(-1*TMath::Cos(fRad2Deg[bin]));
    }

    fFiredPads.resize( fNSector , vector<Int_t>( fNRichSeg , 0 ) );

    return kTRUE;
}



//--------------------------------------------------------------------
// -- Helper functions

void HParticleBtRingF::addressToColRow(const Int_t address, Int_t &sec, Int_t &row, Int_t &col)
{
    //Converts richCal address to sector, column and row numbers

    sec = (Int_t)address/10000;
    row = (Int_t)(address-sec*10000)/100;
    col = address-sec*10000-row*100;
    sec = sec==6 ? 0 : sec;
}

Int_t HParticleBtRingF::correctPhi(const Int_t sec,  const Float_t phi)
{
    //Correct phi angle if phi angle disagrees with sector number
    //Agreement is necessary for transformations from angle to RICH pad plane

    if(sec==0){
        if(phi>120)
            return 120;
        else if(phi<60)
            return 60;
        else
            return phi;
    }else if(sec==1){
        if(phi<120)
            return 120;
        else if(phi>180)
            return 180;
        else
            return phi;
    }else if(sec==2){
        if(phi<180)
            return 180;
        else if(phi>240)
            return 240;
        else
            return phi;
    }else if(sec==3){
        if(phi<240)
            return 240;
        else if(phi>300)
            return 300;
        else
            return phi;
    }else if(sec==4){
        if(phi>360 || (phi>=0 && phi < 60))
            return 360;
        else if(phi<300)
            return 300;
        else
            return phi;
    }else{  //sec==5
        if(phi>300 || phi<0)
            return 0;
        else if(phi >60 && phi<=300)
            return 60;
        else
            return phi;
    }
}

void HParticleBtRingF::sortElements(Double_t &entry1 ,Double_t &entry2)
{
    //Sorts 2 Double_t values by size (ascending)
    if(entry1 >entry2){
        Float_t tmp = entry1;
        entry1      = entry2;
        entry2      = tmp;
    }
}

            
Int_t HParticleBtRingF::getVertexNum(const Float_t vertex)
{
    //Returns vertex number to given z-vertex position

    for(Int_t i=0; i < fNVertex;i++){
        if((vertex>-fVertexPosMin+(i*fVertexStep))&& (vertex<-fVertexPosMax+(i*fVertexStep))) return i;
    }
    if(vertex<=fVertexPosMin)
        return 0;
    else //v>-2.5
        return 14;
}



//--------------------------------------------------------------------
// -- Called once per event

void HParticleBtRingF::fillRichCal(HCategory* catRichCal)
{
    //Fills rich cal hit addresses and charges into vectors
    //Bad hits and hits from bad events are exluded

    Int_t sizeRichCal = catRichCal->getEntries();
    HRichCal* richCal;
    Int_t sec;
    Int_t col;
    Int_t row;
    for(Int_t l=0; l < sizeRichCal;l++){
        richCal = HCategoryManager::getObject(richCal,catRichCal,l);
        if(!richCal)
            continue;
        if(!richCal->getIsCleanedSector() && !richCal->getIsCleanedHigh()){
            addressToColRow(richCal->getAddress(),sec, row, col);
            if(richCal->getCharge() > fChargeLimit[sec]){
                fRichHitAdd.push_back(richCal->getAddress());
                fRichHitCharge.push_back(richCal->getCharge());
            }
             if(row<fRichSegBorderY){
                 fFiredPads[sec][0]++;
             }else{
                 if(col<fRichSegBorderX)
                     fFiredPads[sec][1]++;
                 else
                     fFiredPads[sec][2]++;
             }
        }
    }
}

void HParticleBtRingF::clearData()
{
    //Empties rich and track vectors

    for(UInt_t i = 0; i < fPrediction.size() ; i++ )
        fPrediction[i].clear();
    fPrediction.clear();
    for(UInt_t i = 0; i < fRingMatrix.size() ; i++ )
        fRingMatrix[i].clear();
    fRingMatrix.clear();

    fRichHitAdd.clear();
    fRichHitCharge.clear();
    fTrackTheta.clear();
    fTrackPhi.clear();
    fTrackVertex.clear();
    fTrackSec.clear();
    fTrackPCandIdx.clear();
    fIsGoodTrack.clear();
    fPosXCenter.clear();
    fPosYCenter.clear();

    for(Int_t i = 0; i < fNSector ; i++ ){
        for(Int_t j = 0; j < fNRichSeg ; j++ ){
            fFiredPads[i][j]=0;
        }
    }

}


//--------------------------------------------------------------------
//  -- Called once per track

//Marks pads in area around ring center
void HParticleBtRingF::fillMatrix(Int_t xPad, Int_t yPad, Int_t sec)
{
    //Marks area around ring center
    //Address are stored in ringMatrix vector

    Int_t address    = 0;

    vector <Int_t> ringMatrixVec;
    xPad -= fSizeMatrix*0.5;
    yPad -= fSizeMatrix*0.5;
    address = 10000 * (sec ? sec : 6) + 100 * yPad + xPad;

  
    for(Int_t posX=0; posX <fSizeMatrix; posX++){
        for(Int_t posY=0; posY <fSizeMatrix; posY++){
            if( TMath::Sqrt( TMath::Power(TMath::Abs(fSizeMatrix*0.5-posX),2) + TMath::Power(TMath::Abs(fSizeMatrix*0.5-posY),2) ) < (Float_t)fSizeMatrix/2. ){
                ringMatrixVec.push_back(address+posY*100+posX);
            }
        }
    }
    fRingMatrix.push_back(ringMatrixVec);
}


void HParticleBtRingF::fillPrediction(const HParticleCand* cand, HVertex &vertex, Bool_t isGoodTrack, const Bool_t doAngleCorr)
{
    //Theta and phi angles and vertex position is used to determine track position on pad plane
    //Information about optical ring distortion is used to predict fired pads in defined sigma area
    //Mean values with neighbouring sectors assure a more smooth ring prediction
    //Predicted pads are stored in prediction vector

    HParticleBtAngleTrafo richAngles;

    Int_t address = 0;
    Double_t theta        = cand->getTheta();
    Double_t phi          = cand->getPhi();
    Int_t   sec           = cand->getSector();
    Float_t candZ         = cand->getZ();
    Float_t vertexZ       = vertex.getZ();
    Int_t   vertexNum;
    Double_t radius;
    Double_t sigma;
    Double_t sigmaX0;
    Double_t sigmaY0;
    Double_t radiusSeg0;
    Double_t sigmaSeg0;
    Double_t sigmaX0Seg0;
    Double_t sigmaY0Seg0;
    Double_t radiusNext;
    Double_t sigmaNext;
    Double_t sigmaX0Next;
    Double_t sigmaY0Next;
    Double_t sigNum;
    vector <Int_t> predictionVec;

    phi = correctPhi(sec,phi);

    if(theta>fThetaAngleMin && theta < fThetaAngleMax)
    {
        //Cases if no vertex was reconstructed
        if(vertexZ > -80. && vertexZ < 20.){
            vertexNum     = getVertexNum(vertexZ);
        }else if(candZ > -80. && candZ < 20.){
            vertexZ       = candZ;
            vertexNum     = getVertexNum(candZ);
        }else if(candZ <= -80.){
            vertexZ       = -80.;
            vertexNum     = getVertexNum(-80.);
        }else{
            vertexZ       = 20.;
            vertexNum     = getVertexNum(20.);
        }

     

        //Correct theta and phi angle
        Double_t thetaCor;
        Double_t phiCor;
        if(doAngleCorr){
            if(fAngleCor.alignRichRing(theta,phi, thetaCor, phiCor)){
                theta = thetaCor;
                phi   = phiCor;
            }
        }

        thetaCor = richAngles.zTheta2dTheta(vertexZ,theta, -(phi-fPhiOff[sec]));
        theta-=thetaCor;



        //Transform angles to x,y Padplane coordinates
        Double_t xPad  = richAngles.angles2xPad(theta,-(phi-fPhiOff[sec]));
        Double_t yPad  = richAngles.angles2yPad(theta,-(phi-fPhiOff[sec]));
        Double_t posX  = richAngles.angles2x(theta,-(phi-fPhiOff[sec]));
        Double_t posY  = richAngles.angles2y(theta,-(phi-fPhiOff[sec]));

       
        //z vertex correction
        Double_t yPadCorr = richAngles.zTheta2dYPad(vertexZ,theta,-(phi-fPhiOff2[sec]));
        Double_t xPadCorr = richAngles.zTheta2dXPad();
        Double_t posXCorr = richAngles.zTheta2dY(vertexZ,theta,-(phi-fPhiOff2[sec]));
        Double_t posYCorr = richAngles.zTheta2dX();
        posX += posXCorr;
        posY += posYCorr;
        xPad += xPadCorr;
        yPad += yPadCorr;


        //Coordinates for n sigma areas
        Double_t sigmaX[2]     = {0.,0.};
        Double_t sigmaY[2]     = {0.,0.};
        Double_t padX[2]       = {0.,0.};
        Double_t padY[2]       = {0.,0.};
        Int_t   padXRound[2]   = {0,0};
        Int_t   padYRound[2]   = {0,0};

        sigNum=fNSigma;

        //Loop over different phi angles of the ring
        for( Int_t bin = 0; bin < fNRingSeg; bin++ ) {
          
            //Find neighbouring segments
            Int_t binPrev;
            Int_t binNext;
            if(bin==0){
                binPrev = fNRingSeg-1;
                binNext = 2;
            }
            else if(bin==fNRingSeg-1){
                binPrev = fNRingSeg-2;
                binNext = 0;
            }
            else{
                binPrev = bin-1;
                binNext = bin+1;
            }
            radius  = fPol2DMean[bin][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigma   = fPol2DSigma[bin][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigmaX0 = fRad2DegX[bin]*sigNum*sigma;
            sigmaY0 = fRad2DegY[bin]*sigNum*sigma;


            //Smooth ring prediction area width neighbour sigma and radius
            radiusSeg0  = fPol2DMean[binPrev][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigmaSeg0   = fPol2DSigma[binPrev][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigmaX0Seg0 = fRad2DegX[binPrev]*sigNum*sigmaSeg0;
            sigmaY0Seg0 = fRad2DegY[binPrev]*sigNum*sigmaSeg0;
           
            radiusNext  = fPol2DMean[binNext][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigmaNext   = fPol2DSigma[binNext][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
            sigmaX0Next = fRad2DegX[binNext]*sigNum*sigmaNext;
            sigmaY0Next = fRad2DegY[binNext]*sigNum*sigmaNext;

           
            radius  = (radius + radiusSeg0  + radiusNext  )/3.;
            if(sigma<sigmaSeg0){
                sigma   = (sigma  + sigmaSeg0   + sigmaNext   )/3.;
                sigmaX0 = (sigmaX0+ sigmaX0Seg0 + sigmaX0Next )/3.;
                sigmaY0 = (sigmaY0+ sigmaY0Seg0 + sigmaY0Next )/3.;
            }
       
            //Smooth Overlap region
            if(bin==fNRingSeg-1 || bin==fNRingSeg-2 || bin==fNRingSeg-3){
                radiusSeg0  = fPol2DMean[0][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
                sigmaSeg0   = fPol2DSigma[0][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
                sigmaX0Seg0 = fRad2DegX[0]*sigNum*sigmaSeg0;
                sigmaY0Seg0 = fRad2DegY[0]*sigNum*sigmaSeg0;

                radius  = (radius + radiusSeg0 )/2.;
                if(sigma<sigmaSeg0){
                    sigma   = (sigma  + sigmaSeg0  )/2.;
                    sigmaX0 = (sigmaX0+ sigmaX0Seg0)/2.;
                    sigmaY0 = (sigmaY0+ sigmaY0Seg0)/2.;
                    sigma   = sigmaSeg0  ;
                    sigmaX0 = sigmaX0Seg0;
                    sigmaY0 = sigmaY0Seg0;
                }
            }
            if(bin==0 || bin==1 || bin==2){
                radiusSeg0  = fPol2DMean[fNRingSeg-1][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
                sigmaSeg0   = fPol2DSigma[fNRingSeg-1][vertexNum]->Eval(theta,-(phi-fPhiOff[sec])+fPhiOffsetPar);
                sigmaX0Seg0 = fRad2DegX[fNRingSeg-1]*sigNum*sigmaSeg0;
                sigmaY0Seg0 = fRad2DegY[fNRingSeg-1]*sigNum*sigmaSeg0;

                radius  = (radius + radiusSeg0 )/2.;
                if(sigma<sigmaSeg0){
                    sigma   = (sigma  + sigmaSeg0  )/2.;
                    sigmaX0 = (sigmaX0+ sigmaX0Seg0)/2.;
                    sigmaY0 = (sigmaY0+ sigmaY0Seg0)/2.;
                    sigma   = sigmaSeg0  ;
                    sigmaX0 = sigmaX0Seg0;
                    sigmaY0 = sigmaY0Seg0;
                }
            }


            sigmaX[0] = fRad2DegX[bin]*radius-sigmaX0;
            sigmaY[0] = fRad2DegY[bin]*radius-sigmaY0;
            sigmaX[1] = fRad2DegX[bin]*radius+sigmaX0;
            sigmaY[1] = fRad2DegY[bin]*radius+sigmaY0;
            sortElements(sigmaX[0],sigmaX[1]);
            sortElements(sigmaY[0],sigmaY[1]);
            padX[0]   = richAngles.x2xPad(posX+sigmaX[0]);
            padY[0]   = richAngles.xy2yPad(posX+sigmaX[0],posY+sigmaY[0]);
            padX[1]   = richAngles.x2xPad(posX+sigmaX[1]);
            padY[1]   = richAngles.xy2yPad(posX+sigmaX[1],posY+sigmaY[1]);



            padXRound[0] = TMath::Floor(padX[0]);
            padXRound[1] = TMath::Floor(padX[1]);
            padYRound[0] = TMath::Floor(padY[0]);
            padYRound[1] = TMath::Floor(padY[1]);
        
            //Convert pad position in address and store address in vector
            Int_t padY0=padY[0];
            while(padXRound[0] <= padXRound[1]){
                
                while(padYRound[0] <= padYRound[1]){
                    address = 10000 * (sec ? sec : 6) + 100 * padYRound[0] + padXRound[0];

                    //Fill only if object is not yet stored in vector
                    Bool_t kFound=kFALSE;
                    UInt_t pred=0;
                    while(pred < predictionVec.size() && !kFound){
                        if(predictionVec[pred] == address)
                            kFound=kTRUE;
                        pred++;

                    }
                    if(!kFound)
                        predictionVec.push_back(address);

                    padYRound[0]++;
                }
                padYRound[0]=padY0;
                padXRound[0]++;
            }
        }
        fPrediction.push_back(predictionVec);
        fTrackTheta.push_back(theta);
        fTrackPhi.push_back(phi);
        fTrackVertex.push_back(vertexNum);
        fTrackSec.push_back(sec);
        fTrackPCandIdx.push_back(cand->getIndex());
        fIsGoodTrack.push_back(isGoodTrack);
        fPosXCenter.push_back(xPad);
        fPosYCenter.push_back(yPad);

        fillMatrix(TMath::Floor(xPad),TMath::Floor(yPad),sec);
  
    } else {
        predictionVec.clear();
        predictionVec.push_back(-1);
        fPrediction.push_back(predictionVec);
        fTrackTheta.push_back(-1);
        fTrackPhi.push_back(-1);
        fTrackVertex.push_back(-1);
        fTrackSec.push_back(-1);
        fTrackPCandIdx.push_back(-1);
        fIsGoodTrack.push_back(kFALSE);
        fRingMatrix.push_back(predictionVec);
        fPosXCenter.push_back(-1.);
        fPosYCenter.push_back(-1.);

    }

}



Float_t HParticleBtRingF::getRingMatrix(const Int_t trackNo)
{
    //Returns fraction of  pads fired on ring prediction compared to pads fired in ring matrix around ring center

    Int_t padsFired     = 0;
    Int_t padsFiredPred = 0;
    Float_t chargeFired     = 0;
    Float_t chargeFiredPred = 0;
 
    for(UInt_t i=0; i<fRingMatrix[trackNo].size(); i++){
        for(UInt_t j=0; j<fRichHitAdd.size(); j++){
            if(fRingMatrix[trackNo][i]==fRichHitAdd[j]){
                padsFired++;
                chargeFired +=fRichHitCharge[j];
                for(UInt_t k=0; k<fPrediction[trackNo].size(); k++){
                    if(fRingMatrix[trackNo][i]==fPrediction[trackNo][k]){
                        padsFiredPred++;
                        chargeFiredPred +=fRichHitCharge[j];
                    }
                }
            }
        }
    }
    return (Float_t)chargeFiredPred/(Float_t)chargeFired;
}

Float_t HParticleBtRingF::getTrackTheta(Int_t trackNo)
{
    //Returns theta angle for given trackNo
    if(trackNo<(Int_t)fTrackTheta.size()) {
        return fTrackTheta[trackNo];
    } else return -1;
}

Float_t HParticleBtRingF::getTrackPhi(Int_t trackNo)
{
    //Returns phi angle for given trackNo
    if(trackNo<(Int_t)fTrackPhi.size()) {
        return fTrackPhi[trackNo];
    } else return -1;
}

Int_t HParticleBtRingF::getTrackVertex(Int_t trackNo)
{
    //Returns vertex number for given trackNo
    if(trackNo<(Int_t)fTrackVertex.size()) {
        return fTrackVertex[trackNo];
    } else return -1;
}

Int_t HParticleBtRingF::getTrackSec(Int_t trackNo)
{
    //Returns sector for given trackNo
    if(trackNo<(Int_t)fTrackSec.size()) {
        return fTrackSec[trackNo];
    } else return -1;
}

Float_t HParticleBtRingF::getPosXCenter(Int_t trackNo)
{
    //Returns phi angle for given trackNo
    if(trackNo<(Int_t)fPosXCenter.size()) {
        return fPosXCenter[trackNo];
    } else return -1;
}

Float_t HParticleBtRingF::getPosYCenter(Int_t trackNo)
{
    //Returns phi angle for given trackNo
    if(trackNo<(Int_t)fPosYCenter.size()) {
        return fPosYCenter[trackNo];
    } else return -1;
}

Bool_t HParticleBtRingF::isGoodTrack(Int_t trackNo)
{
    //Returns sector for given trackNo
    if(trackNo<(Int_t)fTrackSec.size()) {
        return fIsGoodTrack[trackNo];
    } else return 0;
}

Bool_t HParticleBtRingF::fillRingInfo(HParticleBtRingInfo* btRingInfo)
{
    //Stores ring,track and rich information in arrays
    //Arrays are stored in HParticleBtRingFInfo for output
    Bool_t kOverflow = kFALSE;

    Int_t   prediction      [128][128]; //tracks,predicted pads
    Int_t   ringMatrix      [128][128];
    Int_t   richHitAdd      [1024];    //fired pads
    Float_t richHitCharge   [1024];

    Float_t trackTheta    [128];
    Float_t trackPhi      [128];
    Int_t trackVertex     [128];
    Int_t trackSec        [128];
    Int_t trackPCandIdx   [128];
    Bool_t  isGoodTrack   [128];
    Float_t posXCenter    [128];
    Float_t posYCenter    [128];


    //Initialization
    for(Int_t i = 0; i < 128; i++ ){
        for(Int_t j = 0; j < 128; j++ ){
            prediction[i][j]  = -1;
            ringMatrix[i][j]  = -1;
        }
    }

    for(Int_t i = 0; i < 1024; i++ ){
        richHitAdd[i]    = -1;
        richHitCharge[i] = -1.;
    }

    for(Int_t i = 0; i < 128; i++ ){
        trackTheta[i]    = -1.;
        trackPhi[i]      = -1.;
        trackVertex[i]   = -1;
        trackSec[i]      = -1;
        trackPCandIdx[i] = -1;
        isGoodTrack[i]   = kFALSE;
        posXCenter[i]    = -1.;
        posYCenter[i]    = -1.;
    }

    //Fill vector to array and check if vector size is too large
    if(fPrediction.size() <= 128){
        for(UInt_t i = 0; i < fPrediction.size(); i++ ){
            if(fPrediction[i].size() <= 128){
                for(UInt_t j = 0; j < fPrediction[i].size(); j++ ){
                    prediction[i][j] = fPrediction[i][j];
                    ringMatrix[i][j] = fRingMatrix[i][j];
                }
            }
            else
                kOverflow = kTRUE;
        }
    }else
        kOverflow = kTRUE;


    if(fRichHitAdd.size() <= 1024){
        for(UInt_t i = 0; i < fRichHitAdd.size(); i++ ){
            richHitAdd[i]     =  fRichHitAdd[i];
            richHitCharge[i]  =  fRichHitCharge[i];
        }
    }else
        kOverflow =kTRUE;
    if(fTrackTheta.size()<=128){
        for(UInt_t i = 0; i < fTrackTheta.size(); i++ ){
            trackTheta[i]    =  fTrackTheta[i];
            trackPhi[i]      =  fTrackPhi[i];
            trackVertex[i]   =  fTrackVertex[i];
            trackSec[i]      =  fTrackSec[i];
            trackPCandIdx[i] =  fTrackPCandIdx[i];
            isGoodTrack[i]   =  fIsGoodTrack[i];
            posXCenter[i]    =  fPosXCenter[i];
            posYCenter[i]    =  fPosYCenter[i];
        }
    }else
        kOverflow =kTRUE;

        btRingInfo->setPrediction       ( prediction   );
        btRingInfo->setRingMatrix       ( ringMatrix   );
        btRingInfo->setRichHitAdd       ( richHitAdd   );
        btRingInfo->setRichHitCharge    ( richHitCharge);
      
        btRingInfo->setTrackTheta       ( trackTheta   );
        btRingInfo->setTrackPhi         ( trackPhi     );
        btRingInfo->setTrackVertex      ( trackVertex  );
        btRingInfo->setTrackSec         ( trackSec     );
        btRingInfo->setTrackPCandIdx    ( trackPCandIdx);
        btRingInfo->setIsGoodTrack      ( isGoodTrack  );
        btRingInfo->setPosXCenter       ( posXCenter   );
        btRingInfo->setPosYCenter       ( posYCenter   );

        return kOverflow;

}

Int_t HParticleBtRingF::plotPrediction(Int_t trackNo = -1)
{
    Int_t returnVal = 0 ;
    for(UInt_t j = 0; j < fPrediction[trackNo].size(); j++ ){
        Bool_t hasHit = kFALSE;
        for(UInt_t k = 0; k <  fRichHitAdd.size(); k++ ){
            if(fRichHitAdd[k] == fPrediction[trackNo][j]){
                hasHit = kTRUE;
                returnVal++;
            }
        }
        cout << fPrediction[trackNo][j] <<  " - " <<"["<< hasHit << "]"  << endl;
    }

return returnVal;
}


void HParticleBtRingF::plotRichHit(Int_t trackNo = -1)
{
    Int_t sec[2];
    Int_t row[2];
    Int_t col[2];
    Int_t plotCounter = 0 ;
    for(UInt_t i = 0; i < fRichHitAdd.size(); i++ ){
        addressToColRow(fRichHitAdd[i]         ,sec[0], row[0], col[0]);
        addressToColRow(fPrediction[trackNo][1],sec[1], row[1], col[1]);
        if(sec[0] == sec[1]){
            cout << fRichHitAdd[i] << " --  " ;
            plotCounter++;
            if(plotCounter%5 == 0)
                cout << endl;
        }
    }
     
}



Bool_t HParticleBtRingF::hasNoisyRichSeg(Bool_t *trackInSec)
{
    Bool_t secUsed[6];
    for(Int_t i=0;i<6;i++) secUsed[0+i] = trackInSec[i];
    Bool_t isNoisy = kFALSE;
    for(Int_t i = 0;  i < fNSector;  i++){
        for(Int_t j = 0;  j < fNRichSeg;  j++){
            if( fFiredPads[i][j] > 100 && secUsed[i])
                isNoisy = kTRUE;
        }
    }
    if(isNoisy)
        cout << "hparticlebtringf: Noisy RICH segment detected and removed" << endl;
    return isNoisy;
}