///////////////////////////////////////////////////////////////////////////////
//  fill RPC branches of GEANT ROOT tree
//
//  Modified on 05/12/2008 by I. Koenig.
//  Modified on 23/12/2009 by D. Gonzalez-Diaz.
//                         Implemented a fast algorithm for filtering
//                         tracks with only valid hits at the RPC box.
//  Modified on 19/07/2012 by Alessio Mangiarotti.
//                         The RPC geometry file in Geant has the lateral
//                         columns exchanged with respect to the accepted
//                         convention and Hydra. As suggested by Jochen
//                         Markert, the columns are swapped here and no
//                         further action is necessary in the digitizer or
//                         in the hit finder.
///////////////////////////////////////////////////////////////////////////////

#include "fillrpc.h"
#include "cfortran.h"
#include "hgeantrpc.h"
#include "hades.h"
#include "hlocation.h"
#include "hcategory.h"
#include "hrecevent.h"
#include "hpartialevent.h"
#include "hgeantkine.h"
#include "hgeantmaxtrk.h"

struct RPC_DEF {    // structure mapped to the common block RPCTUPLE
  Int_t ntrk;
  Int_t rpctrk[MAXTRKRPC];
  Int_t rpcdet[MAXTRKRPC];
  Float_t rpce[MAXTRKRPC];
  Float_t rpcx[MAXTRKRPC];
  Float_t rpcy[MAXTRKRPC];
  Float_t rpcz[MAXTRKRPC];
  Float_t rpctheta[MAXTRKRPC];
  Float_t rpcphi[MAXTRKRPC];
  Float_t rpctof[MAXTRKRPC];
  Float_t rpcmom[MAXTRKRPC];
  Float_t rpclen[MAXTRKRPC];
  Float_t rpcloclen[MAXTRKRPC];
};

#define RPCTUPLE COMMON_BLOCK(RPCTUPLE,rpctuple)
COMMON_BLOCK_DEF(RPC_DEF,RPCTUPLE);
                                                    
void fillrpc(int* sector) {

  HLocation loc;
  HRecEvent* pHGeantEvent = (HRecEvent*)(gHades->getCurrentEvent());
  HPartialEvent* pSimul   = pHGeantEvent->getPartialEvent(catSimul); 
  HCategory* pRpcCat      = pSimul->getCategory(catRpcGeantRaw);         
  HCategory* pKin         = pSimul->getCategory(catGeantKine); 
  HCategory* pKinAux      = pSimul->getCategory(catGeantKine); 

  if(*sector == 1) pRpcCat->Clear();                                

  Int_t index=-1, column=-999, cell=-999, gap=-999;
  Int_t track=0, track_sec=0, track_sec2=0, track_main=0, nbranches=10; //Maximum number of times you
                                                                        //allow to look for father
  Int_t Id, track_tmp;
  Bool_t hasNoDaughterInCell[MAXTRKRPC], hasNot=kFALSE, WasAlreadyTagged=kFALSE;

  Float_t Time_prev=0, Time=0;

  //Beginning of loop for tagging tracks at the entrance of the RPC box
  //with daughters having a geant hit in an RPC cell 

  for(Int_t i=0; i<MAXTRKRPC; i++) hasNoDaughterInCell[i] = kTRUE;

  for(Int_t trk_main=0; trk_main<RPCTUPLE.ntrk; trk_main++) {        

    if(RPCTUPLE.rpcdet[trk_main]<0) continue;     //take only tracks in cell

    track_main = RPCTUPLE.rpctrk[trk_main];    
    Time       = RPCTUPLE.rpctof[trk_main];

    Int_t j          = 0;
    hasNot           = kTRUE;
    WasAlreadyTagged = kFALSE;
    track            = track_main;
    
    while(hasNot && track!=0 && j<nbranches){   
      for(Int_t trk_sec=0; trk_sec<RPCTUPLE.ntrk; trk_sec++) {	

	track_sec = RPCTUPLE.rpctrk[trk_sec];
	Time_prev = RPCTUPLE.rpctof[trk_sec];

	if((track==track_sec)&&!(hasNoDaughterInCell[trk_sec])) {
	  WasAlreadyTagged=kTRUE; 
	  hasNot=kFALSE; 
	  break;
	}
	if((track==track_sec)&&(RPCTUPLE.rpcdet[trk_sec]<0)&&(Time_prev<Time)) {
	  hasNot = kFALSE;	  
	  break;
	}
      }
      track = ((HGeantKine*)pKin->getObject(track-1))->getParentTrack();
      j++;
    }

    if(WasAlreadyTagged) continue;
    
    for(Int_t trk_sec2=0; trk_sec2<RPCTUPLE.ntrk; trk_sec2++) {
      if((RPCTUPLE.rpctrk[trk_sec2]==track_sec) && !hasNot) hasNoDaughterInCell[trk_sec2]   = 0;	  
    }
  }
						      
  //End of loop for tagging

  //Real loop starts for filling starts

  for(Int_t trk=0; trk<RPCTUPLE.ntrk; trk++) {
    Int_t track = RPCTUPLE.rpctrk[trk];
    
    if(hasNoDaughterInCell[trk]) continue;
            
    if(RPCTUPLE.rpcdet[trk]>=0)
      {
	column       =  RPCTUPLE.rpcdet[trk]/(35*4);
	cell         = (RPCTUPLE.rpcdet[trk] - column*(35*4))/4;
	gap          =  RPCTUPLE.rpcdet[trk] - column*(35*4) - cell*4; 
	//
	// ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION
        // Correction for the mistake in the Geant RPC geometry file
	// ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION
	//
        if (column==0) {
	  column=4;
	} else if (column==1) {
	  column=5;
        } else if(column==4) {
	  column=0;
	} else if (column==5) {
          column=1;
	}
	//
        // ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION ATTENTION
	//
    
      }
   
    loc.set(2,(*sector-1),trk); 
           
    HGeantRpc* pRpc = (HGeantRpc*)pRpcCat->getSlot(loc,&index);
    if(pRpc==NULL) {
      printf("\n*** FILLRPC: HGeantRpc location not found ! ***\n");
      printf("sector=%d trk=%d\n",(*sector-1),trk);
      return;
    }

    pRpc = new(pRpc) HGeantRpc;   

    pRpc->setTrack(track);
    pRpc->setHit(RPCTUPLE.rpcx[trk],RPCTUPLE.rpcy[trk],RPCTUPLE.rpcz[trk],
		 RPCTUPLE.rpctof[trk],RPCTUPLE.rpcmom[trk],RPCTUPLE.rpce[trk]);
    pRpc->setIncidence(RPCTUPLE.rpctheta[trk],RPCTUPLE.rpcphi[trk]);
    pRpc->setTLength(RPCTUPLE.rpclen[trk],RPCTUPLE.rpcloclen[trk]);      
    
    if(RPCTUPLE.rpcdet[trk]<0) pRpc->setDetectorID(RPCTUPLE.rpcdet[trk]);
    else pRpc->setAddress(*sector-1, column, cell, gap);

    if(pKin != NULL && track > 0) {  // add TClonesArray index to access from HGeantKine
      ((HGeantKine*)pKin->getObject(track-1))->setRpcHitIndex(index);
      pRpc->sortVariable = RPCTUPLE.rpctof[trk];
    }
  }

  if(pKin != NULL && *sector == 6) { //sort linked lists by sortVariable = RPCTUPLE.rpctof
    Int_t nTracks = pKin->getEntries();
    for(Int_t i=0;i<nTracks;i++) {
      ((HGeantKine*)pKin->getObject(i))->sortRpcHits();
    }
  }
}

// wrapping of the c++ routine as Fortran-subroutine using cfortran.h 
// of Cern-library
FCALLSCSUB1(fillrpc,FILLRPC,fillrpc,PINT)