From d3968b8fa1b1116b5cf883a3214456c41ac91589 Mon Sep 17 00:00:00 2001
From: Florian Uhlig <f.uhlig@gsi.de>
Date: Mon, 6 Jul 2020 12:35:43 +0200
Subject: [PATCH] Add changes in Tof code from dev branch

The merge request ports the commits 16288, 16332, 16335-16336, 16411-16413 and 16509-16511.
---
 reco/detectors/tof/CbmTofCalibrator.cxx       | 129 +++++----
 reco/detectors/tof/CbmTofCalibrator.h         |  18 +-
 reco/detectors/tof/CbmTofEventClusterizer.cxx | 254 +++++++++---------
 reco/detectors/tof/CbmTofFindTracks.cxx       |  66 +++--
 reco/detectors/tof/CbmTofTrackletTools.cxx    |  12 +-
 sim/detectors/tof/CbmTofDigitize.cxx          |  67 +++--
 6 files changed, 309 insertions(+), 237 deletions(-)

diff --git a/reco/detectors/tof/CbmTofCalibrator.cxx b/reco/detectors/tof/CbmTofCalibrator.cxx
index f166cb81..99f179b6 100644
--- a/reco/detectors/tof/CbmTofCalibrator.cxx
+++ b/reco/detectors/tof/CbmTofCalibrator.cxx
@@ -1,7 +1,7 @@
 /** @file CbmTofCalibrator.cxx
  ** @author nh
  ** @date 28.02.2020
- **
+ ** 
  **/
 
 // CBMroot classes and includes
@@ -170,6 +170,7 @@ Bool_t CbmTofCalibrator::CreateCalHist( ){
 				   2*fDigiBdfPar->GetNbChan(iSmType,iRpcId),0,2*fDigiBdfPar->GetNbChan(iSmType,iRpcId),
 				   100, 0., TotMax);
     
+    TSumMax=1.;
     fhCalWalk[iDetIndx].resize( fDigiBdfPar->GetNbChan(iSmType,iRpcId) );
     for( Int_t iCh=0; iCh<fDigiBdfPar->GetNbChan(iSmType,iRpcId); iCh++){
       fhCalWalk[iDetIndx][iCh].resize( 2 );
@@ -187,9 +188,11 @@ Bool_t CbmTofCalibrator::CreateCalHist( ){
 
 void CbmTofCalibrator::FillCalHist( CbmTofTracklet *pTrk){
   // fill deviation histograms on walk level
-  if(pTrk->GetTt() < 0) return;                                   // take tracks with positive velocity only
+  if(pTrk->GetTt() < 0) return;                    // take tracks with positive velocity only
   if( ! pTrk->ContainsAddr( 0x00005006 ) ) return; // request beam counter hit for calibration
-  
+  if (fdR0Lim>0.)  // consider only tracks originating from nominal interaction point
+       if( pTrk->GetR0() > fdR0Lim) return;  
+
   for (Int_t iHit=0; iHit<pTrk->GetNofHits(); iHit++) {
     CbmTofHit* pHit=pTrk->GetTofHitPointer(iHit);
     Int_t iDetId = (pHit->GetAddress() & DetMask);
@@ -212,7 +215,8 @@ void CbmTofCalibrator::FillCalHist( CbmTofTracklet *pTrk){
     Double_t hitpos[3];
     hitpos[0]=pHit->GetX();    hitpos[1]=pHit->GetY();    hitpos[2]=pHit->GetZ();
     Double_t hlocal_p[3];
-//    TGeoNode* cNode= gGeoManager->GetCurrentNode();
+    //TGeoNode* cNode=
+    gGeoManager->GetCurrentNode();
     gGeoManager->MasterToLocal(hitpos, hlocal_p);
     hitpos[0]=pTrk->GetFitX(pHit->GetZ());    hitpos[1]=pTrk->GetFitY(pHit->GetZ());
     Double_t hlocal_f[3];
@@ -236,48 +240,64 @@ void CbmTofCalibrator::FillCalHist( CbmTofTracklet *pTrk){
       
       const CbmTofDigi* tDigi0 = fDigiMan->Get<CbmTofDigi>(iDigInd0);
       Int_t iCh0    = tDigi0->GetChannel();
-      Int_t iSide0 = tDigi0->GetSide();
-      fhCalWalk[iDetIndx][iCh0][iSide0]->Fill(tDigi0->GetTot(),tDigi0->GetTime()-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
-    			                                                 -(1.-2.*tDigi0->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
-					                                 +fTofFindTracks->GetTOff(iDetId)
-					                                );
-
+      Int_t iSide0  = tDigi0->GetSide();
+      LOG(debug)<<"Fill Walk for "<<iDetIndx<<", TSRCS "<<iSmType<<iSm<<iRpc<<iCh0<<iSide0<<", "<<tDigi0<<", "<<pTrk;
+      if( iDetIndx > (Int_t)fhCalWalk.size() ) { LOG(error) << "Invalid DetIndx " << iDetIndx; continue;}
+      if( iCh0 > (Int_t) fhCalWalk[iDetIndx].size() ) { LOG(error) << "Invalid Ch0 " << iCh0; continue;}
+      if( iSide0 > (Int_t) fhCalWalk[iDetIndx][iCh0].size() ) { LOG(error) << "Invalid Side0 " << iSide0; continue;}
+  
+      fhCalWalk[iDetIndx][iCh0][iSide0]->Fill(tDigi0->GetTot(),tDigi0->GetTime()
+                                        +(1.-2.*tDigi0->GetSide())*hlocal_p[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
+                                        -pTrk->GetFitT(pHit->GetZ()) //-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
+                                        +fTofFindTracks->GetTOff(iDetId)
+                                        +2.*(1.-2.*tDigi0->GetSide())*(hlocal_p[1]-hlocal_f[1])/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
+                                        );
       /*      
-      LOG(info)<<"TSRCS "<<iSmType<<iSm<<iRpc<<iCh<<iSide0<<Form(": digo0 %f, ex %f, prop %f, Off %f, res %f",
-								tDigi0->GetTime(),
-								fTrackletTools->GetTexpected(pTrk, iDetId, pHit) ,
-								fTofFindTracks->GetTOff(iDetId),
-								(1.-2.*tDigi0->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc),
-								tDigi0->GetTime()-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
-								-(1.-2.*tDigi0->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc));
-      */					
+      LOG(info)<<"TSRCS "<<iSmType<<iSm<<iRpc<<iCh<<iSide0<<Form(": digi0 %f, ex %f, prop %f, Off %f, res %f",
+                            tDigi0->GetTime(),
+                            fTrackletTools->GetTexpected(pTrk, iDetId, pHit) ,
+                            fTofFindTracks->GetTOff(iDetId),
+                            (1.-2.*tDigi0->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc),
+                            tDigi0->GetTime()-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
+                            -(1.-2.*tDigi0->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc));
+      */
       
       const CbmTofDigi* tDigi1 = fDigiMan->Get<CbmTofDigi>(iDigInd1);
-      Int_t iCh1    = tDigi1->GetChannel();
+      Int_t iCh1   = tDigi1->GetChannel();
       Int_t iSide1 = tDigi1->GetSide();
-      fhCalWalk[iDetIndx][iCh1][iSide1]->Fill(tDigi1->GetTot(),tDigi1->GetTime()-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
-    			                                                 -(1.-2.*tDigi1->GetSide())*hlocal_f[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
-					                                 +fTofFindTracks->GetTOff(iDetId)					      
-					                                );
+      LOG(debug)<<"Fill Walk for "<<iDetIndx<<", TSRCS "<<iSmType<<iSm<<iRpc<<iCh1<<iSide1<<", "<<tDigi1<<", "<<pTrk;
+      if( iCh1 > (Int_t) fhCalWalk[iDetIndx].size() ) { LOG(error) << "Invalid Ch1 " << iCh1; continue;}
+      if( iSide1 > (Int_t) fhCalWalk[iDetIndx][iCh1].size() ) { LOG(error) << "Invalid Side1 " << iSide1; continue;}      
+      fhCalWalk[iDetIndx][iCh1][iSide1]->Fill(tDigi1->GetTot(),tDigi1->GetTime()
+                                        +(1.-2.*tDigi1->GetSide())*hlocal_p[1]/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
+                                        -pTrk->GetFitT(pHit->GetZ()) //-fTrackletTools->GetTexpected(pTrk, iDetId, pHit) 
+                                        +fTofFindTracks->GetTOff(iDetId)
+                                        +2.*(1.-2.*tDigi1->GetSide())*(hlocal_p[1]-hlocal_f[1])/fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc)
+                                        );
     }
-   
   }
 }
 
 Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
   // get current calibration histos
   LOG(info) << "CbmTofCalibrator:: update histos from "
-                 << "file " << CbmTofEventClusterizer::Instance()->GetCalParFileName();
+                   << "file " << CbmTofEventClusterizer::Instance()->GetCalParFileName()
+	           << " with option " << iOpt;
   TFile* fCalParFile = new TFile(CbmTofEventClusterizer::Instance()->GetCalParFileName(),"");
+  if( NULL == fCalParFile ) {
+      LOG(warn) << "Could not open TofClusterizer calibration file, abort Update ";
+      return kFALSE;
+  }
   assert(fCalParFile);
   ReadHist(fCalParFile);
-  
+
+  const Double_t MINCTS=100.; //FIXME, numerical constant in code
   // modify calibration histograms 
   for(Int_t iDetIndx=0; iDetIndx< fDigiBdfPar->GetNbDet(); iDetIndx++){
     Int_t iUniqueId  = fDigiBdfPar->GetDetUId( iDetIndx );
-//    Int_t iSmAddr   = iUniqueId & DetMask; 
+    // Int_t iSmAddr   = iUniqueId & DetMask; 
     Int_t iSmType  = CbmTofAddress::GetSmType( iUniqueId );
-//    Int_t iSm           = CbmTofAddress::GetSmId( iUniqueId );
+    //  Int_t iSm           = CbmTofAddress::GetSmId( iUniqueId );
     Int_t iRpc          = CbmTofAddress::GetRpcId( iUniqueId );
     switch(iOpt) {
     case 0: // none
@@ -285,7 +305,7 @@ Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
     case 1: // update channel mean
       {
 	//LOG(info) << "Update Offsets for TSR "<<iSmType<<iSm<<iRpc;
-//	TProfile *hpP = fhCalPos[iDetIndx]->ProfileX();
+	TProfile *hpP = fhCalPos[iDetIndx]->ProfileX();
 	TProfile *hpT = fhCalTOff[iDetIndx]->ProfileX();
 	TH1* hCalT    = fhCalTOff[iDetIndx]->ProjectionX();
 	//fhCorPos[iDetIndx]->Add((TH1 *)hpP,-1.);
@@ -299,8 +319,12 @@ Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
 	    LOG(info) << Form("Update %s: bin %02d, Cts: %d, Old %f, dev %f, av %f, new %f", fhCorTOff[iDetIndx]->GetName() ,
 			      iBin, (Int_t)dCts, dCorT, dDt, dAvOff, dCorT- dDt + dAvOff); 
 	  }
-	  
-	  fhCorTOff[iDetIndx]->SetBinContent(iBin+1,dCorT+dDt+dAvOff);
+	  Double_t dDp=hpP->GetBinContent(iBin+1);
+	  Double_t dCorP=fhCorPos[iDetIndx]->GetBinContent(iBin+1);
+	  if (dCts > MINCTS) {
+	    fhCorTOff[iDetIndx]->SetBinContent(iBin+1,dCorT+dDt+dAvOff);
+	    fhCorPos[iDetIndx]->SetBinContent(iBin+1,dCorP+dDp);
+	  }
 	}
       }
       break;
@@ -311,8 +335,8 @@ Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
       {
 	for(Int_t iSide=0; iSide<2; iSide++) {
 	  //LOG(info) << "Get walk histo pointer for TSRCS " << iSmType<<iSm<<iRpc<<iCh<<iSide;
-	  TProfile *hpW = fhCalWalk[iDetIndx][iCh][iSide]->ProfileX();        // mean deviation
-	  TH1* hCW        = fhCalWalk[iDetIndx][iCh][iSide]->ProjectionX(); // contributing counts
+      TProfile *hpW = fhCalWalk[iDetIndx][iCh][iSide]->ProfileX();        // mean deviation
+      TH1* hCW      = fhCalWalk[iDetIndx][iCh][iSide]->ProjectionX(); // contributing counts
 
 	  Double_t dCorT=0;
 	  for (Int_t iBin=0; iBin<fhCorTOff[iDetIndx]->GetNbinsX(); iBin++) {
@@ -321,9 +345,27 @@ Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
 	    if (dCts > MinCounts) {
 	      dCorT=hpW->GetBinContent(iBin+1);
 	    }
-	    fhCorWalk[iDetIndx][iCh][iSide]->SetBinContent(iBin+1,dWOff-dCorT); //set new value
+	    fhCorWalk[iDetIndx][iCh][iSide]->SetBinContent(iBin+1,dWOff+dCorT); //set new value
 	  }
-	}
+	  // determine effective/count rate weighted mean
+	  Double_t dMean=0;
+      Double_t dCtsAll=0.;
+	  for (Int_t iBin=0; iBin<fhCorTOff[iDetIndx]->GetNbinsX(); iBin++) {
+	    Double_t dCts=hCW->GetBinContent(iBin+1);
+	    Double_t dWOff= fhCorWalk[iDetIndx][iCh][iSide]->GetBinContent(iBin+1); // current value
+	    if (dCts > MinCounts) {
+          dCtsAll += dCts;
+	      dMean   += dCts*dWOff;
+	    }
+	  }
+	  if(dCtsAll > 0. ) dMean /= dCtsAll;
+      // keep mean value at 0
+      for (Int_t iBin=0; iBin<fhCorTOff[iDetIndx]->GetNbinsX(); iBin++) {
+        Double_t dWOff= fhCorWalk[iDetIndx][iCh][iSide]->GetBinContent(iBin+1); // current value
+        fhCorWalk[iDetIndx][iCh][iSide]->SetBinContent(iBin+1,dWOff-dMean);     //set new value
+      }      
+      
+    }
       }
       break;
     } //switch( iOpt) end  
@@ -331,6 +373,7 @@ Bool_t CbmTofCalibrator::UpdateCalHist( Int_t iOpt){
   
   TFile* fCalParFileNew = new TFile(Form("New_%s",fCalParFile->GetName()),"RECREATE");
   WriteHist(fCalParFileNew);
+  fCalParFileNew->Close();
   
   return kTRUE;
 }
@@ -350,7 +393,7 @@ void CbmTofCalibrator::ReadHist( TFile* fHist){
   
   for(Int_t iDetIndx=0; iDetIndx< fDigiBdfPar->GetNbDet(); iDetIndx++){
     Int_t iUniqueId  = fDigiBdfPar->GetDetUId( iDetIndx );
-//    Int_t iSmAddr   = iUniqueId & DetMask; 
+    //Int_t iSmAddr   = iUniqueId & DetMask; 
     Int_t iSmType  = CbmTofAddress::GetSmType( iUniqueId );
     Int_t iSm           = CbmTofAddress::GetSmId( iUniqueId );
     Int_t iRpc          = CbmTofAddress::GetRpcId( iUniqueId );
@@ -384,29 +427,21 @@ void CbmTofCalibrator::WriteHist( TFile* fHist){
   TDirectory * oldir = gDirectory;
   fHist->cd();
   for(Int_t iDetIndx=0; iDetIndx< fDigiBdfPar->GetNbDet(); iDetIndx++){
-    Int_t iUniqueId  = fDigiBdfPar->GetDetUId( iDetIndx );
-//    Int_t iSmAddr   = iUniqueId & DetMask; 
-    Int_t iSmType  = CbmTofAddress::GetSmType( iUniqueId );
-//    Int_t iSm           = CbmTofAddress::GetSmId( iUniqueId );
-    Int_t iRpc          = CbmTofAddress::GetRpcId( iUniqueId );
-
     fhCorPos[iDetIndx]->Write();
     fhCorTOff[iDetIndx]->Write();
     fhCorTot[iDetIndx]->Write();
     fhCorTotOff[iDetIndx]->Write();
       
-    Int_t iNbCh = fDigiBdfPar->GetNbChan( iSmType, iRpc );
-    fhCorWalk[iDetIndx].resize(iNbCh);
+    Int_t iNbCh = (Int_t)fhCorWalk[iDetIndx].size();
     for( Int_t iCh = 0; iCh < iNbCh; iCh++ )
     {
-      fhCorWalk[iDetIndx][iNbCh].resize(2);
       for(Int_t iSide=0; iSide<2; iSide++){
-	fhCorWalk[iDetIndx][iCh][iSide]->Write();
+        fhCorWalk[iDetIndx][iCh][iSide]->Write();
       }
     }
     
   }
-  gDirectory = oldir;
+  oldir->cd();
 }
 
 ClassImp(CbmTofCalibrator)
diff --git a/reco/detectors/tof/CbmTofCalibrator.h b/reco/detectors/tof/CbmTofCalibrator.h
index 7dd4115f..bc45e5ba 100644
--- a/reco/detectors/tof/CbmTofCalibrator.h
+++ b/reco/detectors/tof/CbmTofCalibrator.h
@@ -54,10 +54,12 @@ class CbmTofCalibrator: public FairTask
 	InitStatus Init();
 	Bool_t  InitParameters();
 	Bool_t  CreateCalHist();
-	void      FillCalHist(CbmTofTracklet *pTrk);
-	Bool_t UpdateCalHist(Int_t iOpt);
-	void ReadHist(TFile *fhFile);		      	
-	void WriteHist(TFile *fhFile);
+	void    FillCalHist(CbmTofTracklet *pTrk);
+	Bool_t  UpdateCalHist(Int_t iOpt);
+	void    ReadHist(TFile *fhFile);		      	
+	void    WriteHist(TFile *fhFile);
+    
+    inline void SetR0Lim(Double_t dVal) {fdR0Lim=dVal;}
 	
  private:
 	CbmDigiManager*               fDigiMan;
@@ -65,9 +67,9 @@ class CbmTofCalibrator: public FairTask
 	CbmTofFindTracks*           fTofFindTracks;
 	CbmTofTrackletTools*     fTrackletTools; 
 
-	CbmTofDigiPar*        fDigiPar;
+	CbmTofDigiPar*    fDigiPar;
 	CbmTofDigiBdfPar* fDigiBdfPar;
-	TClonesArray*           fTofDigiMatchColl;     // TOF Digi Links
+	TClonesArray*     fTofDigiMatchColl;     // TOF Digi Links
 
 	std::vector< TH2* > fhCalPos;    // [nbDet]
 	std::vector< TH2* > fhCalTOff;  // [nbDet]
@@ -82,7 +84,9 @@ class CbmTofCalibrator: public FairTask
 	std::vector< std::vector< std::vector<TH1 *> > >fhCorWalk; // [nbDet][nbCh][nSide]	
 	  
 	std::map<UInt_t, UInt_t> fDetIdIndexMap;
-
+    
+    Double_t fdR0Lim=0.;
+    
 	CbmTofCalibrator(const CbmTofCalibrator&) = delete;
 	CbmTofCalibrator operator=(const CbmTofCalibrator&) = delete;
 	
diff --git a/reco/detectors/tof/CbmTofEventClusterizer.cxx b/reco/detectors/tof/CbmTofEventClusterizer.cxx
index 154a31d6..e531730b 100644
--- a/reco/detectors/tof/CbmTofEventClusterizer.cxx
+++ b/reco/detectors/tof/CbmTofEventClusterizer.cxx
@@ -25,6 +25,7 @@
 #include "CbmEvent.h"
 #include "CbmVertex.h"
 #include "CbmDigiManager.h"
+#include "CbmTofCreateDigiPar.h"      // in tof/TofTools
 
 #include "TTrbHeader.h"
 
@@ -71,11 +72,12 @@ static Int_t    SelMask=DetMask;
 static Double_t fdStartAna10s=0.;
 static Double_t dTLEvt=0.;
 static Int_t iNSpill=0;
+static Int_t iNbTs=0;
 
 const  Double_t fdSpillDuration = 4.;    // in seconds
-const  Double_t fdSpillBreak      = 0.3;  // in seconds
+const  Double_t fdSpillBreak      = 0.9;  // in seconds
 
-static Bool_t bIsMC=kFALSE;
+static Bool_t bAddBeamCounterSideDigi=kTRUE;
 
 //   std::vector< CbmTofPoint* > vPtsRef;
 
@@ -342,25 +344,19 @@ void CbmTofEventClusterizer::SetParContainers()
 
 void CbmTofEventClusterizer::Exec(Option_t* option)
 {
-  /*
-  FairEventHeader* eventHeader = NULL;
-  if ( FairRunAna::Instance() )  eventHeader = FairRunAna::Instance()->GetEventHeader();
-  if( NULL != eventHeader) {
-    LOG(debug) << "MC event entry " <<  eventHeader->GetMCEntryNumber() <<" at t = " << eventHeader->GetEventTime();
-    if(eventHeader->GetMCEntryNumber() > 0) { // first MC event will not be recognized
-      bIsMC=kTRUE; 
-    }
-  }
-  */
+
   if(fTofCalDigiVecOut)  fTofCalDigiVecOut->clear();
   if(fEventsColl) {
-    LOG(debug)<<"CbmTofEventClusterizer::Exec => New time slice with "
-	      << fEventsColl->GetEntriesFast() << " events, "
-	      << fDigiMan->GetNofDigis(ECbmModuleId::kTof) << " TOF digis " ;
+    LOG(info)<<"CbmTofEventClusterizer::Exec => New timeslice " << iNbTs << " with "
+             << fEventsColl->GetEntriesFast() << " events, "
+             << fDigiMan->GetNofDigis(ECbmModuleId::kTof) << " TOF digis " ;
+    iNbTs++;
 
     Int_t iNbHits=0;
     Int_t iNbCalDigis = 0;
     fTofDigiMatchCollOut->Delete(); // costly, FIXME
+    fTofHitsCollOut->Delete();      // costly, FIXME
+    //fTofDigiMatchCollOut->Clear("C"); // not sufficient, memory leak 
     for(Int_t iEvent = 0; iEvent < fEventsColl->GetEntriesFast(); iEvent++)
     {
       CbmEvent* tEvent = dynamic_cast<CbmEvent*>(fEventsColl->At(iEvent));
@@ -374,7 +370,7 @@ void CbmTofEventClusterizer::Exec(Option_t* option)
         fTofDigiVec.push_back(CbmTofDigi(*tDigi));
         //new((*fTofDigisColl)[iNbDigis++]) CbmTofDigi(*tDigi);
       }
-
+    
       ExecEvent(option);
 
       // --- In event-by-event mode: copy caldigis, hits and matches to output array and register them to event
@@ -430,23 +426,6 @@ void CbmTofEventClusterizer::Exec(Option_t* option)
       fTofDigiVec.push_back(CbmTofDigi(*tDigi));
       //new((*fTofDigisColl)[iNbDigis++]) CbmTofDigi(*tDigi);
     }
-
-    if(bIsMC) { // add fake diamond digis, this part should be moved to digitizer
-      //iNbDigis=fTofDigisColl->GetEntriesFast();
-      UInt_t  uChanUId=0x00005006;
-      Double_t dHitTime=gRandom->Gaus(0.,0.04);
-      const Double_t dHitTot=2.;
-      fTofDigiVec.push_back(CbmTofDigi(uChanUId, dHitTime, dHitTot));
-      //CbmTofDigi* tDigi = new CbmTofDigi(uChanUId, dHitTime, dHitTot);
-      //new((*fTofDigisColl)[iNbDigis++]) CbmTofDigi(*tDigi);
- 
-      uChanUId=0x00805006;
-      fTofDigiVec.push_back(CbmTofDigi(uChanUId, dHitTime, dHitTot));
-      //CbmTofDigi* tDigi1 = new CbmTofDigi(uChanUId, dHitTime, dHitTot);
-      //new((*fTofDigisColl)[iNbDigis++]) CbmTofDigi(*tDigi1);
-
-      LOG(debug) << Form("Add fake diamond digis 0x%08x in event mode with t = %7.3f",uChanUId,dHitTime);
-    }
     ExecEvent(option);
   }
 }
@@ -621,12 +600,17 @@ Bool_t   CbmTofEventClusterizer::InitParameters()
    Int_t iGeoVersion = fGeoHandler->Init(isSimulation);
    if( k14a > iGeoVersion )
    {
-      LOG(error)<<"CbmTofEventClusterizer::InitParameters => Only compatible with geometries after v14a !!!"
-                ;
+      LOG(error)<<"CbmTofEventClusterizer::InitParameters => Only compatible with geometries after v14a !!!";
       return kFALSE;
    }
    fTofId = new CbmTofDetectorId_v14a();
    
+   // create digitization parameters from geometry file
+   CbmTofCreateDigiPar* tofDigiPar = new CbmTofCreateDigiPar("TOF Digi Producer",
+                                                        "TOF task");
+   LOG(info) << "Create DigiPar "; 
+   tofDigiPar->Init();
+   
    fDigiPar = (CbmTofDigiPar*) (rtdb->getContainer("CbmTofDigiPar"));
    if( 0 == fDigiPar )
    {
@@ -938,7 +922,7 @@ Bool_t   CbmTofEventClusterizer::LoadGeometry()
      Int_t smodule = fGeoHandler->GetSModule(cellId);
      Int_t module  = fGeoHandler->GetCounter(cellId);
      Int_t cell    = fGeoHandler->GetCell(cellId);
-
+ 
      Double_t x = fChannelInfo->GetX();
      Double_t y = fChannelInfo->GetY();
      Double_t z = fChannelInfo->GetZ();
@@ -957,8 +941,9 @@ Bool_t   CbmTofEventClusterizer::LoadGeometry()
                  ;
       if(icell==0) {
         TGeoHMatrix* cMatrix = gGeoManager->GetCurrentMatrix();
-	fNode->Print();
-	cMatrix->Print();
+	    fNode->Print();
+        fDigiPar->GetNode(cellId)->Print();
+	    cMatrix->Print();
       }
    }
 
@@ -1252,6 +1237,7 @@ Bool_t   CbmTofEventClusterizer::CreateHistos()
                   << ", RpcId "<<iRpcId<<" => UniqueId "<<Form("(0x%08x, 0x%08x)",iUniqueId,iUCellId)
                  <<", dx "<<fChannelInfo->GetSizex()
                  <<", dy "<<fChannelInfo->GetSizey()
+                 <<", z "<<fChannelInfo->GetZ()
                  <<Form(" ChPoi: %p ",fChannelInfo)
                   <<", nbCh "<<fDigiBdfPar->GetNbChan( iSmType, 0 )
                  ;
@@ -1302,7 +1288,7 @@ Bool_t   CbmTofEventClusterizer::CreateHistos()
        fhRpcCluRate10s[iDetIndx] =  new TH1D(
           Form("cl_SmT%01d_sm%03d_rpc%03d_rate10s", iSmType, iSmId, iRpcId ),
           Form("            Clu rate of Rpc #%03d in Sm %03d of type %d in last 10s; Time (s); Rate (Hz)", iRpcId, iSmId, iSmType ),
-	      10000,0.,10.); 
+	      100000,0.,10.); 
 
        fhRpcDTLastHits[iDetIndx] =  new TH1F(
           Form("cl_SmT%01d_sm%03d_rpc%03d_DTLastHits", iSmType, iSmId, iRpcId ),
@@ -1784,7 +1770,8 @@ Bool_t   CbmTofEventClusterizer::FillHistos()
      }
 
      if(fdStartAna10s>0.)
-       fhRpcCluRate10s[iDetIndx]->Fill(dTimeAna10s/1.E9,1./fhRpcCluRate10s[iDetIndx]->GetBinWidth(1));       
+       fhRpcCluRate10s[iDetIndx]->Fill(dTimeAna10s/1.E9,1.);       
+//       fhRpcCluRate10s[iDetIndx]->Fill(dTimeAna10s/1.E9,1./fhRpcCluRate10s[iDetIndx]->GetBinWidth(1));       
      
      if(fdMemoryTime>0. && fvLastHits[iSmType][iSm][iRpc][iCh].size()==0)
        LOG(fatal)<<Form(" <E> hit not stored in memory for TSRC %d%d%d%d",
@@ -2338,7 +2325,6 @@ Bool_t   CbmTofEventClusterizer::FillHistos()
 	     if(pHit->GetZ() < pTrig[iSel]->GetZ()) dZsign[iSel]=-1.;           }
 	   //// look for geometrical match  with selector hit
            if(  iSmType==fiBeamRefType      // to get entries in diamond/BeamRef histos  
-	   //if(  iSmType == 5                  // FIXME, to get entries in diamond histos  
              || TMath::Sqrt(TMath::Power(pHit->GetX()-dzscal*pTrig[iSel]->GetX(),2.)
                            +TMath::Power(pHit->GetY()-dzscal*pTrig[iSel]->GetY(),2.))<fdCaldXdYMax)
            {
@@ -2629,7 +2615,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
       iNent = (Int_t) fhRpcCluAvWalk[iDetIndx]->GetEntries();
      }
      if(0==iNent){
-       LOG(info)<<"WriteHistos: No entries in Walk histos for " 
+       LOG(debug)<<"WriteHistos: No entries in Walk histos for " 
                  <<"Smtype"<<iSmType<<", Sm "<<iSm<<", Rpc "<<iRpc 
                  ;
        // continue;
@@ -2727,7 +2713,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
      }
 
      if(NULL == htempPos_pfx) {
-       LOG(info)<<"WriteHistos: Projections not available, continue " ;
+       LOG(debug)<<"WriteHistos: Projections not available, continue " ;
        continue;
      }
 
@@ -2865,7 +2851,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
           }
           Int_t iNEntries=h2tmp0->GetEntries();
           if(iCh==0)  // condition to print message only once
-          LOG(info)<<Form(" Update Walk correction for SmT %d, Sm %d, Rpc %d, Ch %d, Sel%d: Entries %d",
+          LOG(debug)<<Form(" Update Walk correction for SmT %d, Sm %d, Rpc %d, Ch %d, Sel%d: Entries %d",
                           iSmType,iSm,iRpc,iCh,fCalSel,iNEntries)
                      ;
 
@@ -2892,7 +2878,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
                Double_t dWcor=(((TProfile *)htmp0)->GetBinContent(iWx+1) + ((TProfile *)htmp1)->GetBinContent(iWx+1))*0.5;
                fvCPWalk[iSmType][iSm*iNbRpc+iRpc][iCh][0][iWx]+=dWcor-dWMean;
                fvCPWalk[iSmType][iSm*iNbRpc+iRpc][iCh][1][iWx]+=dWcor-dWMean;
-	       LOG(info) << Form("Walk for TSR %d%d%d%d Tot %d set to %f",iSmType,iSm,iRpc,iCh,iWx, fvCPWalk[iSmType][iSm*iNbRpc+iRpc][iCh][0][iWx]) 
+	       LOG(debug) << Form("Walk for TSR %d%d%d%d Tot %d set to %f",iSmType,iSm,iRpc,iCh,iWx, fvCPWalk[iSmType][iSm*iNbRpc+iRpc][iCh][0][iWx]) 
 			  ;
              }
              break;
@@ -3159,7 +3145,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
         Int_t iNbRpc = fDigiBdfPar->GetNbRpc(  iSmType);
         Int_t iNbCh  = fDigiBdfPar->GetNbChan( iSmType, iRpc );
         if((fCalSmAddr < 0) || (fCalSmAddr != iSmAddr) ){     // select detectors for updating offsets
-         LOG(info)<<"WriteHistos (calMode==3): update Offsets and Gains, keep Walk and DelTof for "
+         LOG(debug)<<"WriteHistos (calMode==3): update Offsets and Gains, keep Walk and DelTof for "
                   <<"Smtype"<<iSmType<<", Sm "<<iSm<<", Rpc "<<iRpc<<" with " <<  iNbCh << " channels "
 		  <<" using selector "<<fCalSel;
          /*
@@ -3308,7 +3294,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
 	      else           TWMean=0.;
 	    }   
 	    if(  htempTOff_px->GetBinContent(iCh+1) > 0. ) 
-	    LOG(info) <<Form("CalibA %d,%2d,%2d: TSRC %d%d%d%d, hits %6.0f, TM %8.3f , TAV  %8.3f, TWM %8.3f, TOff %8.3f,  TOffnew  %8.3f, ",
+	    LOG(debug) <<Form("CalibA %d,%2d,%2d: TSRC %d%d%d%d, hits %6.0f, TM %8.3f , TAV  %8.3f, TWM %8.3f, TOff %8.3f,  TOffnew  %8.3f, ",
 			     fCalMode,fCalSel,fTRefMode,iSmType,iSm,iRpc,iCh,htempTOff_px->GetBinContent(iCh+1), TMean,
 			     ((TProfile *)hAvTOff_pfx)->GetBinContent(iSm*iNbRpc+iRpc+1), TWMean, fvCPTOff[iSmType][iSm*iNbRpc+iRpc][iCh][0],
 			     fvCPTOff[iSmType][iSm*iNbRpc+iRpc][iCh][0]+TMean-TWMean);
@@ -3319,7 +3305,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
           if(htempTOff_px->GetBinContent(iCh+1)>WalkNHmin){
             fvCPTOff[iSmType][iSm*iNbRpc+iRpc][iCh][0] += -dTYOff + TMean;
             fvCPTOff[iSmType][iSm*iNbRpc+iRpc][iCh][1] += +dTYOff + TMean;
-	    LOG(info)<<Form("CalibB %d,%2d,%2d: TSRC %d%d%d%d, hits %6.0f, dTY  %8.3f, TM %8.3f -> new Off %8.3f,%8.3f ",
+	    LOG(debug)<<Form("CalibB %d,%2d,%2d: TSRC %d%d%d%d, hits %6.0f, dTY  %8.3f, TM %8.3f -> new Off %8.3f,%8.3f ",
 			    fCalMode,fCalSel,fTRefMode,iSmType,iSm,iRpc,iCh,htempTOff_px->GetBinContent(iCh+1),
 			  dTYOff,TMean,
 			  fvCPTOff[iSmType][iSm*iNbRpc+iRpc][iCh][0],
@@ -3605,7 +3591,7 @@ Bool_t   CbmTofEventClusterizer::WriteHistos()
         Int_t iNbRpc = fDigiBdfPar->GetNbRpc(  iSmType);
         Int_t iNbCh  = fDigiBdfPar->GetNbChan( iSmType, iRpc );
         if((fCalSmAddr < 0) || (fCalSmAddr != iSmAddr) ){     // select detectors for updating offsets
-         LOG(info)<<"WriteHistos (calMode==5): update Offsets and Gains, keep Walk and DelTof for "
+         LOG(debug)<<"WriteHistos (calMode==5): update Offsets and Gains, keep Walk and DelTof for "
                   <<"Smtype"<<iSmType<<", Sm "<<iSm<<", Rpc "<<iRpc<<" with " <<  iNbCh << " channels "
                    <<" using selector "<<fCalSel
                   ;
@@ -3887,15 +3873,16 @@ Bool_t   CbmTofEventClusterizer::BuildClusters()
     LOG(warning) << "Too many digis in event " << fiNevtBuild  ;
     return kFALSE;
   }
-  if( kFALSE )
+  if( bAddBeamCounterSideDigi )
   {
     // Duplicate type "5" - digis
-    //Int_t iNbDigi=iNbTofDigi; (VF) not used
+    // Int_t iNbDigi=iNbTofDigi; 
     for( Int_t iDigInd = 0; iDigInd < iNbTofDigi; iDigInd++ )
     {
       CbmTofDigi* pDigi = &(fTofDigiVec.at(iDigInd));
       //CbmTofDigi *pDigi = (CbmTofDigi*) fTofDigisColl->At( iDigInd );
       if( pDigi->GetType() == 5 ) {
+        if( pDigi->GetSide() == 1 ) bAddBeamCounterSideDigi=kFALSE; // disable for current data set
         fTofDigiVec.push_back(CbmTofDigi(*pDigi));
         CbmTofDigi* pDigiN = &(fTofDigiVec.back());
         //	 CbmTofDigi *pDigiN  = new((*fTofDigisColl)[iNbDigi++]) CbmTofDigi( *pDigi );
@@ -4434,8 +4421,7 @@ void CbmTofEventClusterizer::CleanLHMemory()
       }
     }
   }
-  LOG(info) << Form("LH cleaning done after %8.0f events",fdEvent)
-	     ;       
+  LOG(info) << Form("LH cleaning done after %8.0f events",fdEvent);       
 }
 
 Bool_t CbmTofEventClusterizer::AddNextChan(Int_t iSmType, Int_t iSm, Int_t iRpc, Int_t iLastChan, Double_t dLastPosX, Double_t dLastPosY, Double_t dLastTime, Double_t dLastTotS){
@@ -4447,6 +4433,8 @@ Bool_t CbmTofEventClusterizer::AddNextChan(Int_t iSmType, Int_t iSm, Int_t iRpc,
   Int_t iDetIndx = fDetIdIndexMap[iDetId];  // Detector Index
 
   Int_t iCh=iLastChan+1;
+  Int_t iChId = CbmTofAddress::GetUniqueAddress(iSm,iRpc,iCh,0,iSmType);
+  
   while( fvDeadStrips[iDetIndx] & ( 1 << iCh )) {
     LOG(debug) << "Skip channel "<<iCh<<" of detector "<<Form("0x%08x",iDetId);
     iCh++; iLastChan++;
@@ -4482,7 +4470,7 @@ Bool_t CbmTofEventClusterizer::AddNextChan(Int_t iSmType, Int_t iSm, Int_t iRpc,
         Double_t dTime = 0.5 * ( xDigiA->GetTime() + xDigiB->GetTime() ) ; 
 	if(TMath::Abs(dTime-dLastTime)<fdMaxTimeDist){
 	  CbmTofDetectorInfo xDetInfo(ECbmModuleId::kTof, iSmType, iSm, iRpc, 0, iCh);
-	  Int_t iChId = fTofId->SetDetectorInfo( xDetInfo );
+	  iChId = fTofId->SetDetectorInfo( xDetInfo );
 	  fChannelInfo = fDigiPar->GetCell( iChId );
 	  gGeoManager->FindNode(fChannelInfo->GetX(),fChannelInfo->GetY(),fChannelInfo->GetZ());
 
@@ -4545,10 +4533,12 @@ Bool_t CbmTofEventClusterizer::AddNextChan(Int_t iSmType, Int_t iSm, Int_t iRpc,
     hitpos_local[1] = (gRandom->Rndm()-0.5)*fChannelInfo->GetSizey()*0.5;
   }
   */
-  Double_t hitpos[3];
-  /*TGeoNode*    cNode   = */gGeoManager->GetCurrentNode();
-  /*TGeoHMatrix* cMatrix = */gGeoManager->GetCurrentMatrix();
-  gGeoManager->LocalToMaster(hitpos_local, hitpos);
+  Double_t hitpos[3]={3*0.};
+  if( 5 != iSmType ) {  // Diamond beam counter always at (0,0,0)
+    /*TGeoNode*    cNode   = */gGeoManager->GetCurrentNode();
+    /*TGeoHMatrix* cMatrix = */gGeoManager->GetCurrentMatrix();
+    gGeoManager->LocalToMaster(hitpos_local, hitpos);
+  }
   TVector3 hitPos(hitpos[0],hitpos[1],hitpos[2]);
   TVector3 hitPosErr(0.5,0.5,0.5);  // FIXME including positioning uncertainty
   Int_t iChm=floor(dLastPosX/fChannelInfo->GetSizex())+iNbCh/2;
@@ -4745,8 +4735,7 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 				  << ", "<<Form("%f",(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][1])->GetTime())
 				  <<", DeltaT " <<(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][1])->GetTime() - 
 			                          (fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0])->GetTime()
-				  <<", array size: " <<  fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size() 
-				  ;
+				  <<", array size: " <<  fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size();
 		       if ( fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][2]->GetSide() 
 			    == fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0]->GetSide() ) {
 			 LOG(debug) << "3 consecutive SameSide Digis! on TSRC "
@@ -4755,8 +4744,7 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 				    << ", "<<(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][1])->GetTime()
 				    <<", DeltaT " <<(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][1])->GetTime() - 
 			                            (fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0])->GetTime()
-				    <<", array size: " <<  fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size() 
-				    ;
+				    <<", array size: " <<  fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size();
 			 fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
 			 fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
 		       }else {
@@ -4771,15 +4759,13 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 			 else {
 			   LOG(debug) 
 			     << Form("Ev %8.0f, digis not properly time ordered, TSRCS %d%d%d%d%d ",
-				     fdEvent,iSmType,iSm,iRpc,iCh,(Int_t)fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0]->GetSide())
-			      ;
+				     fdEvent,iSmType,iSm,iRpc,iCh,(Int_t)fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0]->GetSide());
 			   fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1);
 			   fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1);
 			 }
 		       }
 		     }else{
-		       LOG(debug2)<<"SameSide Erase fStor entries(d) "<<iSmType<<", SR "<<iSm*iNbRpc+iRpc<<", Ch"<<iCh
-				  ;
+		       LOG(debug2)<<"SameSide Erase fStor entries(d) "<<iSmType<<", SR "<<iSm*iNbRpc+iRpc<<", Ch"<<iCh;
 		       fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
 		       fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
 		     }
@@ -4789,12 +4775,10 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 
 		   LOG(debug2) << "digis processing for " 
 			       << Form(" SmT %3d Sm %3d Rpc %3d Ch %3d # %3lu ",iSmType,iSm,iRpc,iCh,
-				       fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size())
-			       ;
+				       fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size());
 		   if(2 > fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size()) {
 		     LOG(debug)<<Form("Leaving digi processing for TSRC %d%d%d%d, size  %3lu",
-				      iSmType,iSm,iRpc,iCh,fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size())
-			       ;
+				      iSmType,iSm,iRpc,iCh,fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size());
 		     break;
 		   }
 		   /* Int_t iLastChId = iChId; // Save Last hit channel*/
@@ -4804,23 +4788,20 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		   iChId = fTofId->SetDetectorInfo( xDetInfo );
 		   Int_t iUCellId=CbmTofAddress::GetUniqueAddress(iSm,iRpc,iCh,0,iSmType);
 		   LOG(debug1)<< Form(" TSRC %d%d%d%d size %3lu ",
-				      iSmType,iSm,iRpc,iCh,fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size())
-			      << Form(" ChId: 0x%08x 0x%08x ",iChId,iUCellId)
-			      ;
+                         iSmType,iSm,iRpc,iCh,fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size())
+                      << Form(" ChId: 0x%08x 0x%08x ",iChId,iUCellId);
 		   fChannelInfo = fDigiPar->GetCell( iChId );
 
 		   if(NULL == fChannelInfo){
 		     LOG(error)<<"CbmTofEventClusterizer::BuildClusters: no geometry info! "
-			       << Form(" %3d %3d %3d %3d 0x%08x 0x%08x ",iSmType, iSm, iRpc, iCh, iChId,iUCellId)
-			       ;
+			       << Form(" %3d %3d %3d %3d 0x%08x 0x%08x ",iSmType, iSm, iRpc, iCh, iChId,iUCellId);
 		     break;
 		   }
 
 		   TGeoNode *fNode=        // prepare local->global trafo
 		     gGeoManager->FindNode(fChannelInfo->GetX(),fChannelInfo->GetY(),fChannelInfo->GetZ());
 		   LOG(debug2)<<Form(" Node at (%6.1f,%6.1f,%6.1f) : %p",
-				     fChannelInfo->GetX(),fChannelInfo->GetY(),fChannelInfo->GetZ(),fNode)
-			      ;
+				     fChannelInfo->GetX(),fChannelInfo->GetY(),fChannelInfo->GetZ(),fNode);
 		   //          fNode->Print();      
 
 		   CbmTofDigi * xDigiA = fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][0];
@@ -4835,8 +4816,7 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		     LOG(debug)<<"CbmTofEventClusterizer::BuildClusters: Diamond hit in "
 			       << iSm <<" with inconsistent digits " 
 			       <<  xDigiA->GetTime() << ", " << xDigiB->GetTime()
-			       << " -> "<<dTimeDif
-			       ;
+			       << " -> "<<dTimeDif;
 		     LOG(debug) << "    "<<xDigiA->ToString();
 		     LOG(debug) << "    "<<xDigiB->ToString();
 		   }
@@ -4847,10 +4827,9 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		     // 0 is the bottom side, 1 is the top side
 		     dPosY = -fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc) * dTimeDif * 0.5;
 
-		   if(TMath::Abs(dPosY) > fChannelInfo->GetSizey() && fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size()>2) {
+		   while(TMath::Abs(dPosY) > fChannelInfo->GetSizey()*fPosYMaxScal && fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size()>2) {
 		     LOG(debug)<<"Hit candidate outside correlation window, check for better possible digis, "
-			       <<" mul "<< fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size()
-			       ;
+			       <<" mul "<< fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].size();
 
 		     CbmTofDigi * xDigiC = fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh][2];
 		     Double_t dPosYN=0.;
@@ -4864,31 +4843,34 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		       dPosYN = -fDigiBdfPar->GetSigVel(iSmType,iSm,iRpc) * dTimeDifN * 0.5;
 
 		     if(TMath::Abs(dPosYN)<TMath::Abs(dPosY)){
-		       LOG(debug)<<"Replace digi on side "<<xDigiC->GetSide()
-				 <<", yPosNext "<<dPosYN
-				 <<" old: "<<dPosY
-				 ;
+		       LOG(debug)<<"Replace digi on side "<<xDigiC->GetSide() <<", yPosNext "<<dPosYN <<" old: "<<dPosY;
 		       dTimeDif=dTimeDifN;
 		       dPosY=dPosYN;
 		       if( xDigiC->GetSide()==xDigiA->GetSide() ) {
-			 xDigiA=xDigiC;
-			 fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
-			 fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
+                 xDigiA=xDigiC;
+                 fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
+                 fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
 		       }else{
-			 xDigiB=xDigiC;
-			 fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(++(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1));
-			 fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(++(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1));
+                 xDigiB=xDigiC;
+                 fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(++(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1));
+                 fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(++(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin()+1));
 		       }
 		     }
-					  
-		   }
-
-		   if(xDigiA->GetSide() == xDigiB->GetSide()){
-		     LOG(fatal)<<"Wrong combinations of digis "
-			       << fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh][0]<<","
-			       << fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh][1]
-			       ;
-		   }
+		     else 
+                 break;
+		   } //while loop end 
+
+           if(xDigiA->GetSide() == xDigiB->GetSide()){
+             LOG(fatal)<<"Wrong combinations of digis "
+                 << fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh][0]<<","
+                 << fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh][1];
+           }
+           
+		   if(TMath::Abs(dPosY) > fChannelInfo->GetSizey()*fPosYMaxScal) {  // remove both digis 
+             fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiExp[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
+             fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].erase(fStorDigiInd[iSmType][iSm*iNbRpc+iRpc][iCh].begin());
+             continue;
+           }
 		   // The "Strip" time is the mean time between each end
 		   dTime    =0.5 * ( xDigiA->GetTime() + xDigiB->GetTime() ) ; 
 
@@ -4963,18 +4945,19 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 			 hitpos_local[1] = (gRandom->Rndm()-0.5)*fChannelInfo->GetSizey();
 		       }
 		       */
-		       Double_t hitpos[3];
-		       TGeoNode*         cNode   = gGeoManager->GetCurrentNode();
-		       /*TGeoHMatrix* cMatrix =*/ gGeoManager->GetCurrentMatrix();
-		       //cNode->Print();
-		       //cMatrix->Print();
-
-		       gGeoManager->LocalToMaster(hitpos_local, hitpos);
-		       LOG(debug1)<<
-			 Form(" LocalToMaster for node %p: (%6.1f,%6.1f,%6.1f) ->(%6.1f,%6.1f,%6.1f)", 
-			      cNode, hitpos_local[0], hitpos_local[1], hitpos_local[2], 
-			      hitpos[0], hitpos[1], hitpos[2])
-				  ;
+		       Double_t hitpos[3]={3*0.};
+               if( 5 != iSmType ) { 
+		         /*TGeoNode*    cNode   =*/ gGeoManager->GetCurrentNode();
+		         /*TGeoHMatrix* cMatrix =*/ gGeoManager->GetCurrentMatrix();
+		         //cNode->Print();
+		         //cMatrix->Print();
+
+		         gGeoManager->LocalToMaster(hitpos_local, hitpos);
+               }
+               LOG(debug1)<<
+                  Form(" LocalToMaster: (%6.1f,%6.1f,%6.1f) ->(%6.1f,%6.1f,%6.1f)", 
+                  hitpos_local[0], hitpos_local[1], hitpos_local[2], 
+                  hitpos[0], hitpos[1], hitpos[2]);
 
 		       TVector3 hitPos(hitpos[0],hitpos[1],hitpos[2]);
 
@@ -5019,15 +5002,13 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 
 		       if(        vDigiIndRef.size() < 2 ){
 			 LOG(warning)<<"Digi refs for Hit "
-				     << fiNbHits<<":        vDigiIndRef.size()"
-				     ;
+				     << fiNbHits<<":        vDigiIndRef.size()";
 		       }                            
 		       if(fiNbHits>0){
 			 CbmTofHit *pHitL = (CbmTofHit*) fTofHitsColl->At(fiNbHits-1);
 			 if(iDetId == pHitL->GetAddress() && dWeightedTime==pHitL->GetTime()){
 			   LOG(debug)<<"Store Hit twice? "
-				     <<" fiNbHits "<<fiNbHits<<", "<<Form("0x%08x",iDetId)
-				     ;
+				     <<" fiNbHits "<<fiNbHits<<", "<<Form("0x%08x",iDetId);
 
 			   for (UInt_t i=0; i<vDigiIndRef.size();i++){
 			     CbmTofDigi *pDigiC = &(fTofCalDigiVec->at(vDigiIndRef.at(i)));
@@ -5037,6 +5018,13 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 			   for (Int_t i=0; i<digiMatchL->GetNofLinks();i++){
 			     CbmLink L0 = digiMatchL->GetLink(i);  
 			     Int_t iDigIndL=L0.GetIndex();
+                 if (iDigIndL >= (Int_t)vDigiIndRef.size() ) {
+                   if( iDetId != fiBeamRefAddr ) {
+                     LOG(warn) << Form("Invalid DigiRefInd for det 0x%08x",iDetId); 
+                     continue; 
+                   }
+                 }
+                 if (vDigiIndRef.at(iDigIndL) >= (Int_t)fTofCalDigiVec->size() ) { LOG(warn) << "Invalid CalDigiInd"; continue; }
 			     CbmTofDigi *pDigiC = &(fTofCalDigiVec->at(vDigiIndRef.at(iDigIndL)));
 			     LOG(debug)<<" DigiL "<<pDigiC->ToString();
 			   }
@@ -5120,8 +5108,8 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		     } // else of if current Digis compatible with last fired chan
 		   } // if( 0 < iNbChanInHit)
 		   else {
-		     LOG(debug)<<Form("1.Hit on TSRC %d%d%d%d, time: %f, PosY %f",iSmType,iSm,iRpc,iCh,dTime,dPosY) 
-			       ;
+		     LOG(debug)<<Form("1.Hit on TSRC %d%d%d%d, time: %f, PosY %f, Tdif %f ",
+                              iSmType,iSm,iRpc,iCh,dTime,dPosY,dTimeDif);
 		     
 		     // first fired strip in this RPC
 		     dWeightedTime = dTime*dTotS;
@@ -5208,18 +5196,18 @@ Bool_t   CbmTofEventClusterizer::BuildHits()
 		 hitpos_local[1] = (gRandom->Rndm()-0.5)*fChannelInfo->GetSizey();
 	       }
 	       */	       
-	       Double_t hitpos[3];
-	       TGeoNode*        cNode= gGeoManager->GetCurrentNode();
-	       /*TGeoHMatrix* cMatrix =*/ gGeoManager->GetCurrentMatrix();
-	       //cNode->Print();
-	       //cMatrix->Print();
-	       
-	       gGeoManager->LocalToMaster(hitpos_local, hitpos);
+	       Double_t hitpos[3]={3*0.};
+           if( 5 != iSmType ) { 
+	         /*TGeoNode*       cNode=*/ gGeoManager->GetCurrentNode();
+	         /*TGeoHMatrix* cMatrix =*/ gGeoManager->GetCurrentMatrix();
+	         //cNode->Print();
+	         //cMatrix->Print();
+	         gGeoManager->LocalToMaster(hitpos_local, hitpos);
+           }
 	       LOG(debug1)<<
-		 Form(" LocalToMaster for V-node %p: (%6.1f,%6.1f,%6.1f) ->(%6.1f,%6.1f,%6.1f)", 
-		      cNode, hitpos_local[0], hitpos_local[1], hitpos_local[2], 
-		      hitpos[0], hitpos[1], hitpos[2])
-			  ;
+		      Form(" LocalToMaster for V-node: (%6.1f,%6.1f,%6.1f) ->(%6.1f,%6.1f,%6.1f)", 
+		      hitpos_local[0], hitpos_local[1], hitpos_local[2], 
+		      hitpos[0], hitpos[1], hitpos[2]);
 	       
 	       TVector3 hitPos(hitpos[0],hitpos[1],hitpos[2]);
 	       // Event errors, not properly done at all for now
@@ -5349,8 +5337,7 @@ Bool_t   CbmTofEventClusterizer::CalibRawDigis()
 	       <<Form("%2d",(Int_t)pDigi->GetChannel())<<" "
 	       <<pDigi->GetSide()<<" "
 	       <<Form("%f",pDigi->GetTime())<<" "
-	       <<pDigi->GetTot()
-	       ;
+	       <<pDigi->GetTot();
 		   
     if(pDigi->GetType()==5 || pDigi->GetType() == 8)   // for Pad counters generate fake digi to mockup a strip
       if(pDigi->GetSide()==1) continue;                // skip one side to avoid double entries
@@ -5386,8 +5373,7 @@ Bool_t   CbmTofEventClusterizer::CalibRawDigis()
 	       <<Form("%2d",(Int_t)pDigi->GetChannel())<<" "
 	       <<pDigi->GetSide()<<" "
 	       <<Form("%f",pDigi->GetTime())<<" "
-	       <<pDigi->GetTot()
-	       ;
+	       <<pDigi->GetTot();
 		   
     if(fbPs2Ns) {
       pCalDigi->SetTime(pCalDigi->GetTime()/1000.);        // for backward compatibility
diff --git a/reco/detectors/tof/CbmTofFindTracks.cxx b/reco/detectors/tof/CbmTofFindTracks.cxx
index bd09ebea..33460e58 100644
--- a/reco/detectors/tof/CbmTofFindTracks.cxx
+++ b/reco/detectors/tof/CbmTofFindTracks.cxx
@@ -18,6 +18,7 @@
 #include "CbmTofDetectorId_v12b.h" // in cbmdata/tof
 #include "CbmTofDetectorId_v14a.h" // in cbmdata/tof
 #include "CbmTofCell.h"       // in tof/TofData
+#include "CbmTofCreateDigiPar.h"      // in tof/TofTools
 #include "CbmTofDigiPar.h"    // in tof/TofParam
 #include "CbmTofDigiBdfPar.h" // in tof/TofParam
 #include "CbmMatch.h"
@@ -285,6 +286,7 @@ InitStatus CbmTofFindTracks::Init()
   // if (fMinNofHits < 1) fMinNofHits=1;
 
     //fill RpcId - map
+    Bool_t bBeamCounter=kFALSE;
     Int_t iRpc=0;
     for (Int_t iCell=0; iCell < fDigiPar->GetNrOfModules(); iCell++){
       Int_t iCellId = fDigiPar->GetCellId(iCell);
@@ -294,12 +296,12 @@ InitStatus CbmTofFindTracks::Init()
 			iRpc, iCellId,
 			fTofId->GetSMType(iCellId), 
                         fTofId->GetSModule(iCellId), 
-                        fTofId->GetCounter(iCellId) )
-		 ;
+                        fTofId->GetCounter(iCellId) );
+        if(fTofId->GetSMType(iCellId) == 5) bBeamCounter=kTRUE;
         fMapRpcIdParInd[iCellId]=iRpc;
-	fRpcAddr.resize(fRpcAddr.size()+1);
-	fRpcAddr.push_back(iCellId);
-	iRpc++;
+	    fRpcAddr.resize(fRpcAddr.size()+1);
+	    fRpcAddr.push_back(iCellId);
+	    iRpc++;
       }
     }
     fStationHMul.resize(fNTofStations+1);
@@ -312,6 +314,10 @@ InitStatus CbmTofFindTracks::Init()
   if(fiCalOpt > 0) {	
    fTofCalibrator  = new CbmTofCalibrator();
    if ( fTofCalibrator->Init() != kSUCCESS ) return kFATAL; 
+   if(bBeamCounter) {
+       fTofCalibrator->SetR0Lim(2.);  // FIXME, hardwired parameter for debugging
+       LOG(info) << "Set CbmTofCalibrator::R0Lim to 2.";
+   }
   }
    
   LOG(info)<<Form("BeamCounter to be used in tracking: 0x%08x",fiBeamCounter);
@@ -428,13 +434,11 @@ Bool_t   CbmTofFindTracks::LoadCalParameter()
 	Int_t iUniqueId = it->first; 
 	CbmTofCell* fChannelInfo   = fDigiPar->GetCell(iUniqueId);
 	if(NULL != fChannelInfo) {
-          Double_t dVal=1.; // FIXME numeric constant in code, default for cosmic 
-	  //if (fiBeamCounter !=-1) 
-	  dVal = fChannelInfo->GetZ() * fTtTarg ; //  use calibration target value
+      Double_t dVal=0.; // FIXME numeric constant in code, default for cosmic 
+	  if (fiBeamCounter != iUniqueId ) 
+        dVal = fChannelInfo->GetZ() * fTtTarg ; //  use calibration target value
 	  fhPullT_Smt_Off->SetBinContent(iDet+1,dVal);
-	  LOG(info)<<Form("Initialize det 0x%08x at %d with TOff %6.2f",
-			  iUniqueId,iDet+1,dVal)
-		   ;
+	  LOG(info)<<Form("Initialize det 0x%08x at %d, z=%f with TOff %6.2f",iUniqueId,iDet+1,fChannelInfo->GetZ(),dVal);
 	}
       }
     }
@@ -509,7 +513,7 @@ Bool_t   CbmTofFindTracks::InitParameters()
       return kFALSE;
    }
 
-   LOG(info)<<"CbmTofFindTTracks::InitParameters: GeoVersion "<<iGeoVersion;
+   LOG(info)<<"CbmTofFindTracks::InitParameters: GeoVersion "<<iGeoVersion;
 
    switch(iGeoVersion){
        case k12b:
@@ -521,6 +525,12 @@ Bool_t   CbmTofFindTracks::InitParameters()
        default:
 	 LOG(error)<<"CbmTofFindTracks::InitParameters: Invalid Detector ID "<<iGeoVersion;
    }
+   
+   // create digitization parameters from geometry file 
+   CbmTofCreateDigiPar* tofDigiPar = new CbmTofCreateDigiPar("TOF Digi Producer","TOF task");
+   LOG(info) << "Create DigiPar ";
+   tofDigiPar->Init();
+   
    return kTRUE;
 }
 // -----  SetParContainers -------------------------------------------------
@@ -707,8 +717,8 @@ Bool_t CbmTofFindTracks::WriteHistos()
 	   if (dFMeanError < 0.05) { // FIXME: hardwired constant 
 	     if(dRMS<RMSmin) dRMS=RMSmin;
 	     if(dRMS>fSIGT*3.0) dRMS=fSIGT*3.;
-	     //if( fRpcAddr[ix] != fiBeamCounter )  // don't correct beam counter position
-	     fhPullT_Smt_Off->SetBinContent(ix+1,dVal);
+	     if( fRpcAddr[ix] != fiBeamCounter )  // don't correct beam counter position
+            fhPullT_Smt_Off->SetBinContent(ix+1,dVal);
 	     fhPullT_Smt_Width->SetBinContent(ix+1,dRMS);
 	   }
 	 }else{
@@ -1182,7 +1192,7 @@ void CbmTofFindTracks::CreateHistograms(){
 
   fhTrklChi2 =  new TH2F(  Form("hTrklChi2"),
 			  Form("Tracklet Chi;  HMul_{Tracklet}; #chi"),
-			  fNTofStations-1, 2, fNTofStations+1, 100, 0, ((CbmTofTrackFinderNN *)fFinder)->GetChiMaxAccept());  
+			  8, 2, 10, 100, 0, ((CbmTofTrackFinderNN *)fFinder)->GetChiMaxAccept());  
   
   fhTrackingTimeNhits  =  new TH2F(  Form("hTrackingTimeNhits"),
 			       Form("Tracking Time; NHits; #Deltat (s)"),
@@ -1190,7 +1200,7 @@ void CbmTofFindTracks::CreateHistograms(){
 
   fhTrklMulNhits =  new TH2F(  Form("hTrklMulNhits"),
 			       Form("Tracklet Multiplicity; NHits; NTracklet"),
-			       100, 0, 200, 20, 0, 20);
+			       150, 0, 150, 25, 0, 25);
 
   fhTrklMulMaxMM =  new TH2F(  Form("hTrklMulMaxMax-1"),
 			       Form("Tracklet Multiplicity; TMulMax; TMulMax-1"),
@@ -1201,31 +1211,31 @@ void CbmTofFindTracks::CreateHistograms(){
 
   fhTrklHMul =  new TH2F(  Form("hTrklHMul"),
 			   Form("Tracklet Hit - Multiplicity; HMul_{Tracklet}; Mul_{HMul}"),
-			   fNTofStations-1, 2, fNTofStations+1, 20, 0, 20);  
+			   8, 2, 10, 20, 0, 20);  
   fhTrklZ0xHMul =  new TH2F(  Form("hTrklZ0xHMul"),
 			   Form("Tracklet Z0x vs. Hit - Multiplicity; HMul_{Tracklet}; Z0x"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -500, 500);
+			   8, 2, 10, 100, -500, 500);
   fhTrklZ0yHMul =  new TH2F(  Form("hTrklZ0yHMul"),
 			   Form("Tracklet Z0y vs. Hit - Multiplicity; HMul_{Tracklet}; Z0y"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -300, 300);
+			   8, 2, 10, 100, -300, 300);
 
   fhTrklTxHMul =  new TH2F(  Form("hTrklTxHMul"),
 			   Form("Tracklet Tx vs. Hit - Multiplicity; HMul_{Tracklet}; Tx"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -0.65, 0.65);
+			   8, 2, 10, 100, -0.65, 0.65);
 
   fhTrklTyHMul =  new TH2F(  Form("hTrklTyHMul"),
 			   Form("Tracklet Ty vs. Hit - Multiplicity; HMul_{Tracklet}; Ty"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -0.65, 0.65);
+			   8, 2, 10, 100, -0.65, 0.65);
   Double_t TTMAX=0.2;
   fhTrklTtHMul =  new TH2F(  Form("hTrklTtHMul"),
 			   Form("Tracklet Tt vs. Hit - Multiplicity; HMul_{Tracklet}; Tt"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -TTMAX, TTMAX);
+			   8, 2, 10, 100, -TTMAX, TTMAX);
   fhTrklVelHMul =  new TH2F(  Form("hTrklVelHMul"),
 			   Form("Tracklet Vel vs. Hit - Multiplicity; HMul_{Tracklet}; v (cm/ns)"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, 0., 50.);
+			   8, 2, 10, 100, 0., 50.);
   fhTrklT0HMul =  new TH2F(  Form("hTrklT0HMul"),
 			   Form("Tracklet T0 vs. Hit - Multiplicity; HMul_{Tracklet}; T0"),
-			   fNTofStations-1, 2, fNTofStations+1, 100, -0.5, 0.5);
+			   8, 2, 10, 100, -0.5, 0.5);
 
   fhTrklT0Mul =  new TH2F( Form("hTrklT0Mul"),
 			   Form("Tracklet #DeltaT0 vs. Trkl - Multiplicity; Mul_{Tracklet}; #Delta(T0)"),
@@ -1535,15 +1545,15 @@ void CbmTofFindTracks::FillHistograms(){
 	vhPullTB[iSt]->Fill(dDTB);
 	vhResidualTBWalk[iSt]->Fill(dTOT,dDTB);
 	vhResidualYWalk[iSt]->Fill(dTOT,dDY);
-	/*
+	
 	fhPullT_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]],dDT);  
 	fhPullX_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]],dDX);  
 	fhPullY_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]],dDY);
-	*/
+	/*
 	fhPullT_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]], fTrackletTools->GetTdif(pTrk,fMapStationRpcId[iSt], pHit) );  
 	fhPullX_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]], fTrackletTools->GetXdif(pTrk,fMapStationRpcId[iSt], pHit) );  
 	fhPullY_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]], fTrackletTools->GetYdif(pTrk,fMapStationRpcId[iSt], pHit) );
-	
+	*/
 	fhPullZ_Smt->Fill((Double_t)fMapRpcIdParInd[fMapStationRpcId[iSt]],dZZ);  
 
 	Double_t dDeltaTt=dTt - fTtTarg;
@@ -1785,7 +1795,7 @@ void CbmTofFindTracks::FillHistograms(){
     }
   }
   if(1)
-    if (fTrackArray->GetEntries() > 20 ) { // temporary  
+    if (fTrackArray->GetEntries() > 25 ) { // temporary  
      LOG(info)<<"Found  high multiplicity of " << fTrackArray->GetEntries() << " in event "<<fiEvent
 	     <<" from "<< fTofHitArray->GetEntries() << " hits " ;
     for (Int_t iTrk=0; iTrk<fTrackArray->GetEntries();iTrk++) {
diff --git a/reco/detectors/tof/CbmTofTrackletTools.cxx b/reco/detectors/tof/CbmTofTrackletTools.cxx
index a3f00a4c..24590eab 100644
--- a/reco/detectors/tof/CbmTofTrackletTools.cxx
+++ b/reco/detectors/tof/CbmTofTrackletTools.cxx
@@ -60,9 +60,9 @@ Double_t* CbmTofTrackletTools::Line3DFit( CbmTofTracklet *pTrk, Int_t iDetId ){
 
 Double_t CbmTofTrackletTools::FitTt( CbmTofTracklet *pTrk, Int_t iDetId ){
   Int_t    nValidHits=0;
-  Int_t    iHit1;
-  Double_t dDist1;
-  //LOG(info) << "FitTt: " << pTrk->GetNofHits() << " hits, exclude " << iDetId;  
+  Int_t    iHit1=0;
+  Double_t dDist1=100.;
+  //LOG(info) << "FitTt: " << pTrk->GetNofHits() << " hits, exclude " << Form("0x%08x",iDetId);  
   Double_t aR[pTrk->GetNofHits()-1];
   Double_t at[pTrk->GetNofHits()-1];
   Double_t ae[pTrk->GetNofHits()-1];
@@ -70,7 +70,7 @@ Double_t CbmTofTrackletTools::FitTt( CbmTofTracklet *pTrk, Int_t iDetId ){
     if ( iDetId == pTrk->GetTofDetIndex(iHit) ) continue;
     if(nValidHits==0){
       iHit1=iHit;
-      //LOG(info) << "Init Dist1";
+      //LOG(info) << "Init Dist1 with " << iHit1;
       Double_t *dPar=Line3DFit(pTrk, iDetId); // spatial fit without iDetId 
       dDist1=pTrk->GetTofHitPointer(iHit1)->GetZ()*TMath::Sqrt(1.+dPar[1]*dPar[1]+dPar[3]*dPar[3]);
       //LOG(info) << "Dist1 = " << dDist1;
@@ -80,8 +80,12 @@ Double_t CbmTofTrackletTools::FitTt( CbmTofTracklet *pTrk, Int_t iDetId ){
     aR[nValidHits]=dDist1+dSign*pTrk->Dist3D(pTrk->GetTofHitPointer(iHit),pTrk->GetTofHitPointer(iHit1));
     at[nValidHits]=pTrk->GetTofHitPointer(iHit)->GetTime();
     ae[nValidHits]=0.1;                         // const timing error, FIXME (?)
+    //LOG(info) << "Init vector index " << nValidHits;
     nValidHits++;    
   }
+  if( nValidHits  < 2 ) return 0.;
+  if( nValidHits == 2 ) return (at[0]-at[1])/(aR[0]-aR[1]);
+  
   Double_t RRsum=0;            //  Values will follow this procedure:
   Double_t Rsum=0;             //  $Rsum=\sum_{i}^{nValidHits}\frac{R_i}{e_i^2}$
   Double_t tsum=0;             //  where e_i will always be the error on the t measurement 
diff --git a/sim/detectors/tof/CbmTofDigitize.cxx b/sim/detectors/tof/CbmTofDigitize.cxx
index 0f030861..953b7086 100644
--- a/sim/detectors/tof/CbmTofDigitize.cxx
+++ b/sim/detectors/tof/CbmTofDigitize.cxx
@@ -46,7 +46,7 @@
 #include "CbmTofDigiPar.h"         // in tof/TofParam
 #include "CbmTofDigiBdfPar.h"      // in tof/TofParam
 #include "CbmTofGeoHandler.h"      // in tof/TofTools
-
+#include "CbmTofCreateDigiPar.h"      // in tof/TofTools
 
 using std::cout;
 using std::endl;
@@ -57,6 +57,7 @@ using std::setw;
 
 // Gauss Integration Constants
 const Int_t    kiNbIntPts = 2;
+Bool_t bFakeBeamCounter=kFALSE;
 
 /************************************************************************************/
 struct CompTimesExp {  
@@ -235,10 +236,8 @@ InitStatus CbmTofDigitize::Init()
   std::cout << std::endl;
   LOG(info) << "=========================================================="
       ;
-  LOG(info) << GetName() << ": Initialisation"  
-      ;
-  if ( fEventMode ) LOG(info) << GetName() << ": Using event mode."
-      ;
+  LOG(info) << GetName() << ": Initialisation";
+  if ( fEventMode ) LOG(info) << GetName() << ": Using event mode.";
 
   // If input file was not set explicitly, use default one
   if ( fsBeamInputFile.IsNull() ) {
@@ -246,7 +245,7 @@ InitStatus CbmTofDigitize::Init()
     fileName += "/parameters/tof/test_bdf_input.root";
     SetInputFileName(fileName);
     LOG(info) << GetName() << ": Using default parameter file "
-        << fileName ;
+              << fileName ;
   }
 
   if( kFALSE == RegisterInputs() )
@@ -409,8 +408,7 @@ Bool_t   CbmTofDigitize::InitParameters()
    Int_t iGeoVersion = fGeoHandler->Init(isSimulation);
    if( k12b > iGeoVersion )
    {
-      LOG(error)<<"CbmTofDigitize::InitParameters => Only compatible with geometries after v12b !!!"
-                ;
+      LOG(error)<<"CbmTofDigitize::InitParameters => Only compatible with geometries after v12b !!!";
       return kFALSE;
    }
 
@@ -426,6 +424,16 @@ Bool_t   CbmTofDigitize::InitParameters()
       default:
          LOG(error)<<"CbmTofDigitize::InitParameters: Invalid Detector ID "<<iGeoVersion;
    }
+   
+   // create digitization parameters from geometry file
+   CbmTofCreateDigiPar* tofDigiPar = new CbmTofCreateDigiPar("TOF Digi Producer",
+                                                        "TOF task");
+   LOG(info) << "Create DigiPar";
+   tofDigiPar->Init();
+   
+   // Printout option for crosscheck
+   LOG(debug)<<"  CbmTofDigitize::LoadBeamtimeValues Digi Par contains "
+             << fDigiPar->GetNrOfModules() << " cells " ;
    return kTRUE;
 }
 Bool_t   CbmTofDigitize::LoadBeamtimeValues()
@@ -434,10 +442,6 @@ Bool_t   CbmTofDigitize::LoadBeamtimeValues()
 
    fDigiBdfPar->SetInputFile(fsBeamInputFile);
 
-   // Printout option for crosscheck
-   LOG(debug)<<"  CbmTofDigitize::LoadBeamtimeValues Digi Par contains "
-             << fDigiPar->GetNrOfModules() << " cells " ;
-
    // Add Param printout only if DEBUG level ON
    if( gLogger->IsLogNeeded( fair::Severity::debug ) )
       fDigiBdfPar->printParams();
@@ -484,6 +488,11 @@ Bool_t   CbmTofDigitize::LoadBeamtimeValues()
       fdChannelGain[iSmType].resize( iNbSm*iNbRpc );
       fvRpcChOffs[iSmType].resize( iNbSm );
 
+      if( iSmType==5 && iNbSm > 0) {
+	bFakeBeamCounter=kTRUE;
+	LOG(info) << "Generate Fake Beam Counter digis";
+      }
+      
       for( Int_t iSm = 0; iSm < iNbSm; iSm++ ){
          fvdSignalVelocityRpc[iSmType][iSm].resize( iNbRpc );
          for( Int_t iRpc = 0; iRpc < iNbRpc; iRpc++ )
@@ -1078,6 +1087,20 @@ Bool_t   CbmTofDigitize::DeleteHistos()
 // TODO: Charge summing up
 Bool_t   CbmTofDigitize::MergeSameChanDigis()
 {
+
+   if( bFakeBeamCounter) {
+     UInt_t  uChanUId=0x00005006;
+     Double_t dHitTime=fCurrentEventTime+gRandom->Gaus(0.,0.04);
+     const Double_t dHitTot=2.;
+     CbmTofDigi *tDigi=new CbmTofDigi(uChanUId, dHitTime, dHitTot);
+     CbmMatch* tMatch = new CbmMatch();
+     SendData(tDigi,tMatch);                        // Send digi to DAQ
+     fiNbDigis++;
+     LOG(debug) << Form("Add fake diamond digis 0x%08x mode with t = %7.3f",uChanUId,dHitTime);
+     //delete tDigi;
+     //delete tMatch;
+   }
+   
    Int_t iNbSmTypes = fDigiBdfPar->GetNbSmTypes();
 
    // loop over each (Smtype, Sm, Rpc, Channel, Side)
@@ -1161,7 +1184,7 @@ Bool_t   CbmTofDigitize::MergeSameChanDigis()
 			      // new digi  will be created and later deleted by CbmDaq.
 			      // The original digi will be deleted below, together with the unused digis from the buffer.
 			      CbmTofDigi* digi = new CbmTofDigi(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].first));
-                  CbmMatch* match = new CbmMatch(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].second));
+			      CbmMatch* match = new CbmMatch(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].second));
 
 			      digi->SetTime(digi->GetTime() * fdDigiTimeConvFactor + fCurrentEventTime);  // ns->ps
 			      SendData(digi, match);                        // Send digi to DAQ
@@ -1200,7 +1223,7 @@ Bool_t   CbmTofDigitize::MergeSameChanDigis()
 			// new digi  will be created and later deleted by CbmDaq.
 			// The original digi will be deleted below, together with the unused digis from the buffer.
 			CbmTofDigi* digi = new CbmTofDigi(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].first));
-            CbmMatch* match = new CbmMatch(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].second));
+			CbmMatch* match = new CbmMatch(*(fStorDigi[iSmType][iSm*iNbRpc + iRpc][iNbSides*iCh+iSide][iDigi0].second));
 			digi->SetTime(digi->GetTime() * fdDigiTimeConvFactor + fCurrentEventTime);  // ns->ps
 			SendData(digi, match);                        // Send digi to DAQ
                         fiNbDigis++;
@@ -1243,7 +1266,7 @@ Bool_t   CbmTofDigitize::MergeSameChanDigis()
                   } // for each (Ch, Side) pair
             } // for each (Sm, rpc) pair
       } // for( Int_t iSmType = 0; iSmType < iNbSmTypes; iSmType++ )
-
+      
    return kTRUE;
 }
 /************************************************************************************/
@@ -1552,9 +1575,15 @@ Bool_t   CbmTofDigitize::DigitizeDirectClusterSize()
 		   ;
       gGeoManager->GetCurrentNode();
       gGeoManager->GetCurrentMatrix();
+      
       Double_t poipos[3]={vPntPos.X(),vPntPos.Y(),vPntPos.Z()};
-      Double_t poipos_local[3];
-      gGeoManager->MasterToLocal(poipos, poipos_local);
+      Double_t poipos_local[3]; 
+      Double_t poipos_local_man[3]; 
+      gGeoManager->MasterToLocal(poipos, poipos_local_man);
+      fDigiPar->GetNode(iChanId)->MasterToLocal(poipos, poipos_local);
+      for (Int_t i=0; i<3; i++) 
+          if ( poipos_local[i] != poipos_local_man[i] ) 
+              LOG(fatal) << "Inconsistent M2L result " << i << ": "<< poipos_local[i] << " != " << poipos_local_man[i];
 
       if( 1 == iChType)
       {
@@ -1778,6 +1807,8 @@ Bool_t   CbmTofDigitize::DigitizeDirectClusterSize()
                            + ( poipos_local[1] )
 #endif
                              /fvdSignalVelocityRpc[iSmType][iSM][iRpc];
+                             LOG(debug)<<"Create DigiA TSRC " << iSmType << iSM << iRpc << iStripInd 
+                                       <<Form("at %f, ypos %f",dTimeA,poipos_local[1]);
                      CbmTofDigi * tofDigi = new CbmTofDigi( iSM, iRpc, iStripInd, dTimeA,
                            dStripCharge*fdChannelGain[iSmType][iSM*iNbRpc + iRpc][2*iStripInd+1]/2.0,
                            1, iSmType );
@@ -1800,6 +1831,8 @@ Bool_t   CbmTofDigitize::DigitizeDirectClusterSize()
                            - ( poipos_local[1] )
 #endif
                              /fvdSignalVelocityRpc[iSmType][iSM][iRpc];
+                             LOG(debug)<<"Create DigiB TSRC " << iSmType << iSM << iRpc << iStripInd 
+                                       <<Form("at %f, ypos %f",dTimeB,poipos_local[1]);
                      CbmTofDigi * tofDigi = new CbmTofDigi( iSM, iRpc, iStripInd, dTimeB,
                            dStripCharge*fdChannelGain[iSmType][iSM*iNbRpc + iRpc][2*iStripInd]/2.0,
                            0, iSmType );
-- 
GitLab