#include "CbmTrdDigitizerPRF.h"

#include "CbmTrdDigiPar.h"
#include "CbmTrdPoint.h"
#include "CbmTrdDigi.h"
#include "CbmTrdDigiMatch.h"
#include "CbmTrdModule.h"
#include "CbmTrdRadiator.h"
#include "CbmTrdGeoHandler.h"

#include "CbmMCTrack.h"

#include "FairRootManager.h"
#include "FairRunAna.h"
#include "FairRuntimeDb.h"

#include "TRandom.h"
#include "TMath.h"
#include "TVector3.h"
#include "TF1.h"
#include "TH2F.h"
#include "TH1F.h"
#include "TProfile.h"
#include "TCanvas.h"
#include "TImage.h"
#include "TClonesArray.h"
#include "TGeoManager.h"
#include "TStopwatch.h"

#include "TPRegexp.h"

//#include "sqlite/sqlite3.h" // used for the lookup table option

#include <iostream>
#include <cmath>
using std::cout;
using std::endl;
using std::pair;
using std::make_pair;
using std::map;
using std::vector;
using std::max;
using std::fabs;
// ---- Default constructor -------------------------------------------
CbmTrdDigitizerPRF::CbmTrdDigitizerPRF()
  : FairTask("TrdCluster"),
    fTime(-1.),
    fModuleType(-1),
    fModuleCopy(-1),
    fModuleID(-1),
    fMCindex(-1),
    fnRow(-1),
    fnCol(-1),
    fnSec(-1),
    fModuleId(-1),
    fLayerId(-1),
    fTrdPoints(NULL),
    fDigis(new TClonesArray("CbmTrdDigi")),
    fDigiMatchs(NULL),
    fMCStacks(NULL),
    fDigiPar(NULL),
    fModuleInfo(NULL),
    fRadiators(NULL),
    fDigiMap(),
    fGeoHandler(new CbmTrdGeoHandler())
{
}
// --------------------------------------------------------------------

// ---- Constructor ----------------------------------------------------
CbmTrdDigitizerPRF::CbmTrdDigitizerPRF(const char *name, const char *title,
					     CbmTrdRadiator *radiator)
  :FairTask(name),
   fTime(-1.),
   fModuleType(-1),
   fModuleCopy(-1),
   fModuleID(-1),
   fMCindex(-1),
    fnRow(-1),
    fnCol(-1),
    fnSec(-1),
    fModuleId(-1),
    fLayerId(-1),
   fTrdPoints(NULL),
   fDigis(NULL),
   fDigiMatchs(NULL),
   fMCStacks(NULL),
   fDigiPar(NULL),
   fModuleInfo(NULL),
   fRadiators(radiator),
   fDigiMap(),
   fGeoHandler(new CbmTrdGeoHandler())
{
}
// --------------------------------------------------------------------

// ---- Destructor ----------------------------------------------------
CbmTrdDigitizerPRF::~CbmTrdDigitizerPRF()
{
  FairRootManager *ioman =FairRootManager::Instance();
  ioman->Write();
  //fDigis->Clear("C");
  fDigis->Delete();
  delete fDigis;
  //fDigiMatchs->Clear("C");
  fDigiMatchs->Delete();
  delete fDigiMatchs;
  fTrdPoints->Clear("C");
  fTrdPoints->Delete();
  delete fTrdPoints;
  fMCStacks->Clear("C");
  fMCStacks->Delete();
  delete fMCStacks;
  if(fDigiPar)
    delete fDigiPar;
  if(fModuleInfo)
    delete fModuleInfo;
  //delete fdigi;
  //delete fdigiMatch;
}
// --------------------------------------------------------------------

// ----  Initialisation  ----------------------------------------------
void CbmTrdDigitizerPRF::SetParContainers()
{
    cout<<" * CbmTrdDigitizerPRF * :: SetParContainers() "<<endl;


    // Get Base Container
    FairRunAna* ana = FairRunAna::Instance();
    FairRuntimeDb* rtdb=ana->GetRuntimeDb();

    fDigiPar = (CbmTrdDigiPar*)
               (rtdb->getContainer("CbmTrdDigiPar"));

}
// --------------------------------------------------------------------

// ---- ReInit  -------------------------------------------------------
InitStatus CbmTrdDigitizerPRF::ReInit(){

  cout<<" * CbmTrdClusterizer * :: ReInit() "<<endl;


  FairRunAna* ana = FairRunAna::Instance();
  FairRuntimeDb* rtdb=ana->GetRuntimeDb();

  fDigiPar = (CbmTrdDigiPar*)
      (rtdb->getContainer("CbmTrdDigiPar"));
  
  return kSUCCESS;
}
// --------------------------------------------------------------------

// ---- Init ----------------------------------------------------------
InitStatus CbmTrdDigitizerPRF::Init()
{

  cout<<" * CbmTrdDigitizerPRF * :: Init() "<<endl;

  FairRootManager *ioman = FairRootManager::Instance();
  if ( ! ioman ) Fatal("Init", "No FairRootManager");
    
  fTrdPoints=(TClonesArray *)ioman->GetObject("TrdPoint"); 
  if ( ! fTrdPoints ) {
    cout << "-W CbmTrdClusterFast::Init: No TrdPoints array!" << endl;
    cout << "                            Task will be inactive" << endl;
    return kERROR;
  }

  fMCStacks = (TClonesArray*)ioman->GetObject("MCTrack");
  if ( ! fMCStacks ) {
    cout << "-W CbmTrdClusterFast::Init: No MCTrack array!" << endl;
    cout << "                            Task will be inactive" << endl;
    return kERROR;
  }

  fDigis = new TClonesArray("CbmTrdDigi", 100);
  ioman->Register("TrdDigi","TRD Digis",fDigis,kTRUE);

  fDigiMatchs = new TClonesArray("CbmTrdDigiMatch", 100);
  ioman->Register("TrdDigiMatch","TRD Digis",fDigiMatchs,kTRUE);

  fGeoHandler->Init();

  fRadiators->Init();

  return kSUCCESS;

}
// --------------------------------------------------------------------


// ---- Exec ----------------------------------------------------------
void CbmTrdDigitizerPRF::Exec(Option_t * option)
{
  Bool_t Debug = false;
  TStopwatch timer;
  timer.Start();
  cout << "================CbmTrdDigitizerPRF=====================" << endl;
 
  CbmTrdPoint *point = NULL;
  CbmMCTrack *track = NULL;

  TVector3 mom;
  Double_t ELoss;
  Double_t ELossTR;     // TR energy loss for e- & e+
  Double_t ELossdEdX;   // ionization energy loss
  Double_t time;
  Int_t nCount = 0;
  Int_t nEntries = fTrdPoints->GetEntries();
  Int_t pointId;
  cout << " Found " << nEntries << " MC-Points in Buffer of TRD" << endl;
  for (Int_t j = 0; j < nEntries ; j++ ){
    pointId = j;
    ELoss = 0.0;
    ELossTR = 0.0;
    ELossdEdX = 0.0;
    //cout << " " << j << endl;
    point = (CbmTrdPoint*) fTrdPoints->At(j);
    if(NULL == point) {
      cout << " no point found " << endl;
    }
    if(NULL == point) continue;

    point->Momentum(mom);
    fMCindex=point->GetTrackID();
    track = (CbmMCTrack*) fMCStacks->At(fMCindex);
    if(NULL == track) {
      cout << " no point found " << endl;
    }

    if(NULL == track) continue;

    // Add TR energy loss for electrons ans positrons
    Int_t pdgCode = track->GetPdgCode();
    if(TMath::Abs(pdgCode) == 11){ 
      ELossTR = fRadiators->GetTR(mom);
    }
    ELossdEdX = point->GetEnergyLoss();
    ELoss = ELossTR + ELossdEdX;

    fTime = point->GetTime();

    Double_t point_in[3];
    Double_t point_out[3];

    point_in[0] = point->GetXIn();
    point_in[1] = point->GetYIn();
    point_in[2] = point->GetZIn();

    point_out[0] = point->GetXOut();
    point_out[1] = point->GetYOut();
    point_out[2] = point->GetZOut();

    // Find node corresponding to the point in the center between entrance and exit MC-point coordinates
    gGeoManager->FindNode((point_out[0] + point_in[0]) / 2, 
			  (point_out[1] + point_in[1]) / 2, 
			  (point_out[2] + point_in[2]) / 2
			  );

 
 
    const Int_t nCluster = 1;//TODO: as function of track length in active volume
   
    const Double_t *global_point = gGeoManager->GetCurrentPoint();
    Double_t local_point[3];

    gGeoManager->MasterToLocal(global_point, local_point);

    fModuleID = point->GetDetectorID();
    fModuleInfo = fDigiPar->GetModule(fModuleID);
    // Form address which contains layerId, moduleId, sectorId, columnId and rowId.
    // Address in the MC point contains only information about layerId and moduleId.
     fLayerId = CbmTrdAddress::GetLayerId(point->GetDetectorID());
     fModuleId = CbmTrdAddress::GetModuleId(point->GetDetectorID());


    SplitTrackPath(point, ELoss);

     /*
     // To be deleted after debuging--------------------------------- start

    Double_t local_point_out[3];
    Double_t local_point_in[3]; // entrace point coordinates in local module cs
    Double_t cluster_pos[3]; // cluster position in local module coordinate system
    Double_t cluster_delta[3]; // vector pointing in MC-track direction with length of one path slice within chamber volume to creat n cluster
    gGeoManager->MasterToLocal(point_in,  local_point_in);
    gGeoManager->MasterToLocal(point_out, local_point_out);
    for (Int_t i = 0; i < 3; i++){
      cluster_delta[i] = (local_point_out[i] - local_point_in[i]) / (nCluster + 1);
      cluster_pos[i] = local_point_in[i] + cluster_delta[i] * nCluster;    
      //printf("(%i):  point_in %f  point_out %f  global_point %f  local_point %f   local_point_in %f  \n",i,point_in[i],point_out[i],global_point[i],local_point[i],local_point_in[i]);
      //printf("   (%i):  %4.f = %4.f + %4.f * %i\n",i,cluster_pos[i],local_point_in[i],cluster_delta[i], nCluster);
    }
    if (Debug)
      if ( fabs(cluster_pos[1]) > fModuleInfo->GetSizeY() )
	{
	  printf("ori %d\n", fModuleInfo->GetOrientation() );
	  printf("%f %f\n", cluster_pos[1], fModuleInfo->GetSizeY() );
	  printf("%f %f\n", local_point_in[1], cluster_delta[1]  );
	  printf("%f %f\n", point_out[1], point_in[1]  );
	}
    fModuleInfo->ProjectPositionToNextAnodeWire(cluster_pos);

    // Form address which contains layerId, moduleId, sectorId, columnId and rowId.
    // Address in the MC point contains only information about layerId and moduleId.
    //Int_t fLayerId = CbmTrdAddress::GetLayerId(point->GetDetectorID());
    //Int_t fModuleId = CbmTrdAddress::GetModuleId(point->GetDetectorID());
   

    const Int_t nCol = fModuleInfo->GetNofColumns();
    const Int_t nRow = fModuleInfo->GetNofRows();

    fnRow = fModuleInfo->GetNofRows();
    fnCol = fModuleInfo->GetNofColumns();
    fnSec = fModuleInfo->GetNofSectors();
    Int_t sectorId(0), columnId(0), rowId(0); // globale row and column ids
    Int_t preSecRows (0);
    fModuleInfo->GetPadInfo(point, sectorId, columnId, rowId);
    if (sectorId < 0) {
      continue;
    }
    if (sectorId > 0){
      for (Int_t iSec = 0; iSec < sectorId; iSec++){
	rowId += fModuleInfo->GetNofRowsInSector(iSec); // sector row id -> global row
	preSecRows += fModuleInfo->GetNofRowsInSector(iSec); // all full sector rows up to the one we are inside
      }
    } 
    Double_t displacement_x(0), displacement_y(0);
    Double_t h = fModuleInfo->GetAnodeWireToPadPlaneDistance();
    Double_t W(fModuleInfo->GetPadSizeX(sectorId)), H(fModuleInfo->GetPadSizeY(sectorId));
    fModuleInfo->TransformToLocalPad(cluster_pos, displacement_x, displacement_y);
    //Define pad-plane area of interest: maximum 7 columns and 3 rows
    const Int_t maxCol(7), maxRow(3);
    Int_t startCol(columnId-maxCol/2), stopCol(columnId+maxCol/2+1), startRow(rowId-maxRow/2), stopRow(rowId+maxRow/2+1), startSec(0);
    if (startCol < 0 )  
      startCol = 0;
    if (stopCol  > nCol) 
      stopCol  = nCol;
    if (startRow < 0)   
      startRow = 0;
    if (stopRow  > nRow) 
      stopRow  = nRow;

    if (startRow > fModuleInfo->GetNofRowsInSector(0))
      startSec = 1;
    if (startRow > fModuleInfo->GetNofRowsInSector(0) + fModuleInfo->GetNofRowsInSector(1))
      startSec = 2;
    Int_t iSec = startSec;
    // Scan pad-plane area to determine charge distribution
    Double_t sum = 0;
    for (Int_t iRow = startRow; iRow < stopRow; iRow++){
      Int_t globalRow = iRow - preSecRows;
      if (globalRow < 0){
	globalRow = iRow;
      }
      if (Debug) printf("[%2i]: ",iRow);
      for (Int_t iCol = startCol; iCol < stopCol; iCol++){
	Int_t address = CbmTrdAddress::GetAddress(fLayerId, fModuleId, sectorId, globalRow, iCol); 
	//printf("[%f, %f, %f]\n",cluster_pos[0],cluster_pos[1],cluster_pos[2]);
	//fModuleInfo->TransformToLocalPad(cluster_pos, displacement_x, displacement_y);
	Double_t chargeFraction = 0;
	if (rowId == iRow && columnId == iCol)
	  chargeFraction = CalcPRF(displacement_x, W, h) * 
	    CalcPRF(displacement_y, H, h);
	else
	  chargeFraction = CalcPRF((iCol - columnId) * W - displacement_x, W, h) * 
	    CalcPRF((iRow - rowId) * H - displacement_y, H, h);
	sum += chargeFraction;
	map<Int_t, pair<CbmTrdDigi*, CbmTrdDigiMatch*> >::iterator it = fDigiMap.find(address);
	if (it == fDigiMap.end()) { // Pixel not yet in map -> Add new pixel
	  fDigiMap[address] = make_pair(new CbmTrdDigi(address, chargeFraction*ELoss, fTime), new CbmTrdDigiMatch(pointId));
	} else { // Pixel already in map -> Add charge
	  CbmTrdDigi* digi = it->second.first;
	  digi->AddCharge(ELoss);
	  digi->SetTime(max(fTime, digi->GetTime()));
	  CbmTrdDigiMatch* digiMatch = it->second.second;
	  digiMatch->AddPoint(pointId);
	}
	if (Debug) {
	  if (rowId == iRow && columnId == iCol)
	    printf("(%3i) ",iCol);
	  else
	    printf(" %3i  ",iCol);
	  if (iCol == nCol-1)
	    cout << "|";
	}
      }
      if (Debug) {
	if (iRow == nRow-1)
	  cout << endl << "-------------------------------------" << endl;
	else
	  cout << endl;
      }
    }
    if (Debug) cout << endl;
    if (sum < 0.95){
      LOG(DEBUG2) << "Summed charge: " << std::setprecision(5) << sum << "  hit:(" << columnId << ", " << rowId <<")   max:(" << nCol-1 << ", " << nRow-1 << ")" << FairLogger::endl;
      //printf("Summed charge: %f  hit:(%i, %i) max:(%i, %i)\n",sum,columnId,rowId,nCol-1,nRow-1);
    }
  }
  // To be deleted after debuging---------------------------------stop
  */
  }
  // Fill data from internally used stl map into output TClonesArray
  Int_t iDigi = 0;
  map<Int_t, pair<CbmTrdDigi*, CbmTrdDigiMatch*> >::iterator it;
  for (it = fDigiMap.begin() ; it != fDigiMap.end(); it++) {
    CbmTrdDigi* digi = it->second.first;
    new ((*fDigis)[iDigi]) CbmTrdDigi(*digi);
    CbmTrdDigiMatch* digiMatch = it->second.second;
    new ((*fDigiMatchs)[iDigi]) CbmTrdDigiMatch(*digiMatch);
    delete digi;
    delete digiMatch;
    iDigi++;
  }
  fDigiMap.clear();
  printf("\n   Created %i TrdDigis\n",iDigi);
  timer.Stop();
  Double_t rtime = timer.RealTime();
  Double_t ctime = timer.CpuTime();

  printf("\n\n******************** Reading Test  **********************\n");
  printf("   RealTime=%f seconds, CpuTime=%f seconds\n",rtime,ctime);
  printf("*********************************************************\n\n");

  //fIntegralTest->Draw();
}
// --------------------------------------------------------------------
Double_t CbmTrdDigitizerPRF::CalcPRF(Double_t x, Double_t W, Double_t h)
{
  Float_t K3 = 0.525; 
  Double_t SqrtK3 = sqrt(K3);

  return fabs(
	      -1. / (2. * atan(SqrtK3)) * (
					   atan(SqrtK3 * tanh(TMath::Pi() * (-2. + SqrtK3 ) * (W + 2.* x) / (8.* h) )) + 
					   atan(SqrtK3 * tanh(TMath::Pi() * (-2. + SqrtK3 ) * (W - 2.* x) / (8.* h) ))
					   )
	      );
}
// --------------------------------------------------------------------
void CbmTrdDigitizerPRF::ScanPadPlane(const Double_t* local_point, Double_t clusterELoss)
{
  // TODO: scan full surface but sum up if pad is outside module untill reenter the active pad-plane-> charge conservation
  Int_t sectorId(-1), columnId(-1), rowId(-1), preSecRows(0);
  fModuleInfo->GetPadInfo( local_point, sectorId, columnId, rowId);// TODO: David add this function to TrdModule;
  if (sectorId < 0 && columnId < 0 && rowId < 0) {
    return;
  } else {
    Double_t displacement_x(0), displacement_y(0);
    Double_t h = fModuleInfo->GetAnodeWireToPadPlaneDistance();
    Double_t W(fModuleInfo->GetPadSizeX(sectorId)), H(fModuleInfo->GetPadSizeY(sectorId));
    fModuleInfo->TransformToLocalPad(local_point, displacement_x, displacement_y);
    displacement_x *= 10.; //cm->mm
    displacement_y *= 10.; //cm->mm
    const Int_t maxCol(7), maxRow(3);
    //Estimate starting column and row and limits due to chamber dimensions
    Int_t startCol(columnId-maxCol/2), stopCol(columnId+maxCol/2+1), startRow(rowId-maxRow/2), stopRow(rowId+maxRow/2+1), startSec(0);
    if (startCol < 0 )  
      startCol = 0;
    if (stopCol  > fnCol) 
      stopCol  = fnCol;
    if (startRow < 0)   
      startRow = 0;
    if (stopRow  > fnRow) 
      stopRow  = fnRow;

    if (startRow > fModuleInfo->GetNofRowsInSector(0))
      startSec = 1;
    if (startRow > fModuleInfo->GetNofRowsInSector(0) + fModuleInfo->GetNofRowsInSector(1))
      startSec = 2;
    Int_t iSec = startSec; // is different from sectorId (id of the maximum pad in the center of the scanned region) if lowest pad row is in a different sector compared to pad max
    // Scan pad-plane area to determine charge distribution
    Double_t sum = 0;
    for (Int_t iRow = startRow; iRow < stopRow; iRow++){
      Int_t globalRow = iRow - preSecRows;
      if (globalRow < 0){
	globalRow = iRow;
      }
      //if (Debug) printf("[%2i]: ",iRow);
      for (Int_t iCol = startCol; iCol < stopCol; iCol++){
	Int_t address = CbmTrdAddress::GetAddress(fLayerId, fModuleId, sectorId, globalRow, iCol); 
	//printf("[%f, %f, %f]\n",cluster_pos[0],cluster_pos[1],cluster_pos[2]);
	//fModuleInfo->TransformToLocalPad(cluster_pos, displacement_x, displacement_y);
	Double_t chargeFraction = 0;
	if (rowId == iRow && columnId == iCol)
	  chargeFraction = CalcPRF(displacement_x, W, h) * 
	    CalcPRF(displacement_y, H, h);
	else
	  chargeFraction = CalcPRF((iCol - columnId) * W - displacement_x, W, h) * 
	    CalcPRF((iRow - rowId) * H - displacement_y, H, h);
	sum += chargeFraction;
	AddDigi(fMCindex, address, Double_t(chargeFraction * clusterELoss), fTime);
	/*if (Debug) {
	  if (rowId == iRow && columnId == iCol)
	    printf("(%3i) ",iCol);
	  else
	    printf(" %3i  ",iCol);
	  if (iCol == fnCol-1)
	    cout << "|";
	    }*/
      }
      /*
      if (Debug) {
	if (iRow == fnRow-1)
	  cout << endl << "-------------------------------------" << endl;
	else
	  cout << endl;
      }
      */
    }
    //if (Debug) cout << endl;
    if (sum < 0.95){
      LOG(DEBUG2) << "Summed charge: " << std::setprecision(5) << sum << "  hit:(" << columnId << ", " << rowId <<")   max:(" << fnCol-1 << ", " << fnRow-1 << ")" << FairLogger::endl;
      //printf("Summed charge: %f  hit:(%i, %i) max:(%i, %i)\n",sum,columnId,rowId,nCol-1,nRow-1);
    }
  }  
}
// --------------------------------------------------------------------
void CbmTrdDigitizerPRF::SplitTrackPath(const CbmTrdPoint* point, Double_t ELoss)
{
  const Double_t nClusterPerCm = 1.0;//TODO: as function of track length in active volume
  Double_t point_in[3] = {
    point->GetXIn(),
    point->GetYIn(),
    point->GetZIn()
};
  Double_t point_out[3] = {
    point->GetXOut(),
    point->GetYOut(),
    point->GetZOut()
  };
  Double_t local_point_out[3];// exit point coordinates in local module cs
  Double_t local_point_in[3]; // entrace point coordinates in local module cs
  gGeoManager->MasterToLocal(point_in,  local_point_in);
  gGeoManager->MasterToLocal(point_out, local_point_out);

  Double_t cluster_pos[3];   // cluster position in local module coordinate system
  Double_t cluster_delta[3]; // vector pointing in MC-track direction with length of one path slice within chamber volume to creat n cluster

  Double_t trackLength = 0;

  for (Int_t i = 0; i < 3; i++){
    cluster_delta[i] = (local_point_out[i] - local_point_in[i]);
    trackLength += cluster_delta[i] * cluster_delta[i];
  }
  //trackLength = TMath::Sqrt(trackLength);
  const Int_t nCluster = TMath::Sqrt(trackLength) / nClusterPerCm + 0.9;// Track length threshold of minimum 0.1cm track length in gas volume
  if (nCluster < 1){
    //LOG(DEBUG2) << "Summed charge: " << std::setprecision(5) << sum << "  hit:(" << columnId << ", " << rowId <<")   max:(" << fnCol-1 << ", " << fnRow-1 << ")" << FairLogger::endl;
    printf("nCluster: %i   track length: %fcm  nCluster/cm: %3.1f  ELoss: %EGeV\n",nCluster, TMath::Sqrt(trackLength), nClusterPerCm, ELoss);
    return;
  }
  Double_t clusterELoss = ELoss / nCluster;
  for (Int_t iCluster = 1; iCluster <= nCluster; iCluster++){
    for (Int_t i = 0; i < 3; i++){
      cluster_pos[i] = local_point_in[i] + cluster_delta[i] * iCluster;  
    }

    fModuleInfo->ProjectPositionToNextAnodeWire(cluster_pos);

    ScanPadPlane(cluster_pos, clusterELoss);
  }
}
// --------------------------------------------------------------------
void CbmTrdDigitizerPRF::AddDigi(Int_t pointId, Int_t address, Double_t charge, Double_t time)
{
  map<Int_t, pair<CbmTrdDigi*, CbmTrdDigiMatch*> >::iterator it = fDigiMap.find(address);
  if (it == fDigiMap.end()) { // Pixel not yet in map -> Add new pixel
    fDigiMap[address] = make_pair(new CbmTrdDigi(address, charge, time), new CbmTrdDigiMatch(pointId));
  } else { // Pixel already in map -> Add charge
    CbmTrdDigi* digi = it->second.first;
    digi->AddCharge(charge);
    digi->SetTime(max(time, digi->GetTime()));
    CbmTrdDigiMatch* digiMatch = it->second.second;
    digiMatch->AddPoint(pointId);
  }
}
// ---- FinishTask-----------------------------------------------------
void CbmTrdDigitizerPRF::FinishEvent()
{  
  fDigis->Delete();
  fDigiMatchs->Delete();
  fTrdPoints->Clear();
  fMCStacks->Clear();
}
// ---- Register ------------------------------------------------------
void CbmTrdDigitizerPRF::Register()
{
  FairRootManager::Instance()->Register("TrdDigi","Trd Digi", fDigis, kTRUE);
  FairRootManager::Instance()->Register("TrdDigiMatch","Trd Digi Match", fDigiMatchs, kTRUE);
}
// --------------------------------------------------------------------

ClassImp(CbmTrdDigitizerPRF)