// -------------------------------------------------------------------------
// -----                PndHypStripHitProducer            -----
// -------------------------------------------------------------------------


#include "TClonesArray.h"
#include "TArrayD.h"
#include "TVector2.h"
#include "TGeoManager.h"

#include "FairRootManager.h"
#include "PndHypStripHitProducer.h"
#include "PndHypHit.h"
#include "PndHypPoint.h"
#include "FairRunAna.h"
#include "FairRuntimeDb.h"
#include "FairGeoNode.h"
//#include "FairRuntimeDb.h"
#include "FairGeoNode.h"
#include "FairGeoVector.h"
//#include "StringVector.h"
#include "PndHypCalcStrip.h"
#include "PndHypDigiStrip.h"
//#include "PndHypDigiStrip.h"
//enum SensorSide { kTOP, kBOTTOM };
// -----   Default constructor   -------------------------------------------
PndHypStripHitProducer::PndHypStripHitProducer() :
FairTask("Hyp Strip Hit Producer")
{
  fBranchName   = "HypPoint";
  
  /*ftopPitch = 0.;
   fbotPitch = 0.;
   forient = 0.;
   fskew = 0.;
   fTopAnchor = TVector2(0,0);
   fBotAnchor = TVector2(0,0);
   fnrTopFE = 0;
   fnrBotFE = 0;
   fnrFECh = 0;
   fthreshold = 0.;
   fnoise = 0.;*/
  
  fOverrideParams = false;
  // stripHits = 0;
  //  fHitArray  = new TClonesArray("PndHypHit");
  //	fStripArray	= new TClonesArray("PndHypStripHit");
}
// -------------------------------------------------------------------------

PndHypStripHitProducer::PndHypStripHitProducer(Double_t topPitch, Double_t botPitch,
                                               Double_t ori, Double_t skew,
                                               TVector2 topAnchor, TVector2 botAnchor,
                                               Int_t nrTopFE, Int_t nrBotFE, Int_t nrFECh,
                                               Double_t threshold, Double_t noise,
                                               TString sensorType, TString feType) :
FairTask("Hyp Strip Hit Producer")
{
  fBranchName   = "HypPoint";
  // ftopPitch = topPitch;
  //   fbotPitch = botPitch;
  //   forient = ori;
  //   fskew = skew;
  //   fTopAnchor = topAnchor;
  //   fBotAnchor = botAnchor;
  //   fnrTopFE = nrTopFE;
  //   fnrBotFE = nrBotFE;
  //   fnrFECh = nrFECh;
  //   fthreshold = threshold;
  //   fnoise = noise;
  
  fOverrideParams = true;
  SetParamSet(topPitch,botPitch,ori,skew,topAnchor,botAnchor,nrTopFE,nrBotFE,
              nrFECh,threshold,noise,"Rect",feType);
  //stripHits = 0;
  std::cout << "Hyp Strip Hit Producer initiated" << std::endl;
}
// -------------------------------------------------------------------------

void PndHypStripHitProducer::SetParamSet(Double_t topPitch, Double_t botPitch,
                                         Double_t ori, Double_t skew,
                                         TVector2 topAnchor, TVector2 botAnchor,
                                         Int_t nrTopFE, Int_t nrBotFE, Int_t nrFECh,
                                         Double_t threshold, Double_t noise,
                                         TString sensorType, TString feType)
{
  FairRunAna* ana = FairRunAna::Instance();
  FairRootManager* ioman = FairRootManager::Instance();
  if ( fDigiPar==0 ) SetParContainers();
  if (fOverrideParams){
    if (sensorType.Contains("Rect"))fCurrentDigiPar = fDigiPar;
    else if (sensorType.Contains("Trap")) fCurrentDigiPar = fDigiPar;
    fCurrentDigiPar->SetTopPitch(topPitch);
    fCurrentDigiPar->SetBotPitch(botPitch);
    fCurrentDigiPar->SetSkew(skew);
    fCurrentDigiPar->SetOrient(ori);
    fCurrentDigiPar->SetTopAnchor(topAnchor);
    fCurrentDigiPar->SetBotAnchor(botAnchor);
    fCurrentDigiPar->SetNrFECh(nrFECh);
    fCurrentDigiPar->SetNrTopFE(nrTopFE);
    fCurrentDigiPar->SetNrBotFE(nrBotFE);
    fCurrentDigiPar->SetThreshold(threshold);
    fCurrentDigiPar->SetNoise(noise);
    fCurrentDigiPar->SetSensType(sensorType);
    fCurrentDigiPar->SetFeType(feType);
    fCurrentDigiPar->setChanged();
    fCurrentDigiPar->setInputVersion(ana->GetRunId(),1);
  }
 	
 	
}

// -----   Destructor   ----------------------------------------------------
PndHypStripHitProducer::~PndHypStripHitProducer()
{ 
  delete fGeoH;
}
// -------------------------------------------------------------------------

// -----   Initialization  of Parameter Containers -------------------------
void PndHypStripHitProducer::SetParContainers()
{
  // Get Base Container
  
  FairRunAna* ana = FairRunAna::Instance();
  FairRuntimeDb* rtdb=ana->GetRuntimeDb();
  //fGeoPar = (PndGeoHypPar*)(rtdb->getContainer("PndGeoHypPar"));
  fDigiPar = (PndHypStripDigiPar*)(rtdb->getContainer("PndHypStripDigiPar"));
  //fGeoMappingPar = (PndHypGeoMappingPar*)(rtdb->getContainer("PndHypGeoMappingPar"));
  std::cout << "fdigi par init" << std::endl;
}

InitStatus PndHypStripHitProducer::ReInit()
{
  /*
   FairRunAna* ana = FairRunAna::Instance();
   FairRuntimeDb* rtdb=ana->GetRuntimeDb();
   //fGeoPar = (PndGeoHypPar*)(rtdb->getContainer("PndGeoHypPar"));
   fDigiPar=(PndHypStripDigiPar*)(rtdb->getContainer("PndHypStripDigiPar"));
   //fGeoMappingPar = (PndHypGeoMappingPar*)(rtdb->getContainer("PndHypGeoMappingPar"));
   std::cout << "fdigi par 2 init" << std::endl;
   */
  SetParContainers();
  return kSUCCESS;
  
}


// -----   Public method Init   --------------------------------------------
InitStatus PndHypStripHitProducer::Init()
{
  FairRunAna* ana = FairRunAna::Instance();
  FairRootManager* ioman = FairRootManager::Instance();
  
  fGeoH = new PndHypGeoHandling(gGeoManager);
  
  if ( ! ioman ) 
  {
    std::cout << "-E- PndHypStripHitProducer::Init: "
    << "RootManager not instantiated!" << std::endl;
    return kFATAL;
  }
  
  fPointArray = (TClonesArray*) ioman->GetObject(fBranchName);
  if ( ! fPointArray )
  {
    std::cout << "-W- PndHypStripHitProducer::Init: "
    << "No HypPoint array!" << std::endl;
    return kERROR;
  }
  
  // Create and register output array
  //fHitArray = new TClonesArray("PndHypHit");
  //ioman->Register("HypHit", "Hyp", fHitArray, kTRUE);
  
  // Create and register output array
  fStripArray = new TClonesArray("PndHypDigiStrip");
  ioman->Register("HypStripDigis", "Hyp", fStripArray, kTRUE);
  
  // Create and register parameter array
  //  fStripArray = new TClonesArray("PndHypDigiPar");
  //  ioman->Register("HypDigiParam", "Hyp", fDigiPar, kTRUE);
  
  //  SetParContainers();
  
  std::cout << "-I- PndHypStripHitProducer: Initialisation successfull" << std::endl;
  //std::cout << "-I- Test: " << fGeoMappingPar->detIdPairArray->GetEntries() << std::endl;
  
  
  if (!fDigiPar) {
    std::cout<<"-E- MvdStripHitProducer: DigiPar Container does not exist!"<<std::endl;
    return kERROR;
  }
  // if (fOverrideParams)
  //        {
  // 	 fDigiPar->SetTopPitch(ftopPitch);
  // 	 fDigiPar->SetBotPitch(fbotPitch);
  // 	 fDigiPar->SetSkew(fskew);
  // 	 fDigiPar->SetOrient(forient);
  // 	 fDigiPar->SetTopAnchor(fTopAnchor);
  // 	 fDigiPar->SetBotAnchor(fBotAnchor);
  // 	 fDigiPar->SetNrFECh(fnrFECh);
  // 	 fDigiPar->SetNrTopFE(fnrTopFE);
  // 	 fDigiPar->SetNrBotFE(fnrBotFE);
  // 	 fDigiPar->SetThreshold(fthreshold);
  // 	 fDigiPar->SetNoise(fnoise);
  // 	 fDigiPar->setChanged();
  // 	 fDigiPar->setInputVersion(ana->GetRunId(),1);
  
  //        }
  
  fDigiPar->print();
  
  fStripCalcTop = new PndHypCalcStrip(fDigiPar, kTOP);
  fStripCalcBot = new PndHypCalcStrip(fDigiPar, kBOTTOM);
  fStripCalcTop->SetVerboseLevel(fVerbose);
  fStripCalcBot->SetVerboseLevel(fVerbose);
  
  return kSUCCESS;
}
// -------------------------------------------------------------------------



// -----   Public method Exec   --------------------------------------------
void PndHypStripHitProducer::Exec(Option_t* opt)
{
  // Reset output array
  //if ( ! fHitArray ) 
  // Fatal("Exec", "No HitArray");
  
  //fHitArray->Delete();
  
  fStripArray->Delete();
  
  //  fFeStripArray->Delete();
  
  // Declare some variables
  PndHypPoint *point = NULL;
  
  Int_t detID = 0,       // Detector ID
  trackID = 0;     // Track index
  
  // Loop over PndHypPoints
  Int_t 
  nPoints = fPointArray->GetEntriesFast();
  if (fVerbose > 0){
    std::cout<<" Nr of Points: "<<nPoints<<std::endl;
  }
  
  
  Int_t iStrip = 0;
  //Int_t iFeStrip = 0;
  Int_t nSiL =-1;
  
  /*PndHypCalcStrip StripCalcTop(fDigiPar->topPitch, fDigiPar->orient,
   fDigiPar->topNrFE * fDigiPar->feChannels , fDigiPar->feChannels,
   fDigiPar->topAnchor,
   fDigiPar->threshold, fDigiPar->noise);
   
   PndHypCalcStrip StripCalcBot(fDigiPar->botPitch, fDigiPar->orient + fDigiPar->skew,
   fDigiPar->botNrFE * fDigiPar->feChannels, fDigiPar->feChannels,
   fDigiPar->botAnchor,
   fDigiPar->threshold, fDigiPar->noise);
   
   PndHypCalcStrip StripCalcTop(ftopPitch, forient,
   fnrTopFE * fnrFECh , fnrFECh,
   fTopAnchor,
   fthreshold, fnoise);
   
   PndHypCalcStrip StripCalcBot(fbotPitch, forient + fskew,
   fnrBotFE * fnrFECh, fnrFECh,
   fBotAnchor,
   fthreshold, fnoise);
   */
  
  for (Int_t iPoint = 0; iPoint < nPoints; iPoint++) 
	{
	  point = (PndHypPoint*) fPointArray->At(iPoint);
	  //if (!((point->GetDetName()).Contains("Strip"))) continue;
	  //if ((point->GetDetName()).Contains("Disk")) continue;
	  //      if ((point->GetDetName()).Contains("SideChips")) continue;
	  
	  //if (fVerbose > 1){
	  //std::cout<<"***** Strip Hit ******"<<std::endl;
	  //std::cout<<" DetName : "<<point->GetDetName()<<std::endl;
	  //}
	  if ( ! point){
	    std::cout<< "No Point!" << std::endl;
	    continue;
	  }
	  if (fVerbose > 1){
	    std::cout << "****Global Point: " << std::endl;
	    point->Print("");
	  }
	  //FairGeoVector posInL, posOutL, meanPos, meanPosL;
	  //GetLocalHitPoints(point, posInL, posOutL);
	  TVector3 posIn,posOut;
	  point->PositionIn(posIn);
	  point->PositionOut(posOut);
	  
	  TVector3 posInL = fGeoH->MasterToLocalId(posIn,point->GetDetName());
	  TVector3 posOutL = fGeoH->MasterToLocalId(posOut,point->GetDetName());
	  
	  if (fVerbose > 1){  
	    std::cout << "Local Hits: " << std::endl;
	    
	  }
    
	  if (fVerbose > 1){
	    
	    std::cout << "Energy: " << point->GetEnergyLoss() << std::endl;
	  }
	  
	  //TVector3 size = GetSensorDimensions(point->GetDetName().Data());
	  TVector3 size = GetSensorDimensions(fGeoH->GetPath(point->GetDetName()).Data());
	  
	  if (fVerbose > 2) std::cout<< "Sensor Size : x="
      <<size.x()*2.<<" , y="
      <<size.y()*2.<<" , z ="
      <<size.z()*2.<<std::endl;
	  
	  TVector2 anchor(-size.x(),size.y());
	  
	  // Top Side
	  if (fVerbose > 0)
	    std::cout  << "Top Side: " << std::endl;
	  // PndHypCalcStrip StripCalc(ftopPitch, forient, fnrTopFE * fnrFECh , 
	  // 				 //FairGeoVector(-size.x(),-size.y(),0.),
	  // 				 TVector2(0.,0.),
	  // 				 fthreshold, fnoise);
	  
	  //StripCalcTop.SetVerboseLevel(fVerbose);
	  Double_t sttop = ((size.x()*2)/fDigiPar->GetTopPitch());
	  
	  Int_t rfetop = Int_t(sttop+0.5);
	  Int_t rtfe = Int_t((rfetop/128)+0.5);
    
	  if(rfetop>=rtfe*128)fStripCalcTop->SetreStrip(rtfe*128);
	  else fStripCalcTop->SetreStrip(rfetop);
	  std::cout<<" tof fe "<<Int_t(rfetop/128)<<" total strips "<<rfetop<<std::endl;
	  
	  nSiL = point->GetVolumeID();
    
	  fStripCalcTop->SetAnchor(anchor);
	  // Calculate a cluster of Strips fired
	  std::vector<PndHypStrip> topStrips =
    fStripCalcTop->GetStrips(posInL.X(), posInL.Y(), posInL.Z(),
                             posOutL.X(), posOutL.Y(), posOutL.Z(),
                             point->GetEnergyLoss(),nSiL);
	  
	  if (topStrips.size() != 0){
	    //stripHits += topStrips.size();
	    if (fVerbose > 0)
	      std::cout  << "SensorStrips: top " << topStrips.size()<<std::endl;
	    
	    trackID = point->GetTrackID();
	    detID   = point->GetVolumeID();
      
	    //      for(uint i = 0; i < myStrips.size(); i++)
	    //iStrip = 0;
	    Int_t sp;SensorSide si;
	    
	    
	    
	    for(std::vector<PndHypStrip>::const_iterator it=topStrips.begin();
          it!= topStrips.end();it++)
        
      {
        
        new ((*fStripArray)[iStrip]) PndHypDigiStrip(iPoint,detID,
                                                     point->GetDetName(),
                                                     fStripCalcTop->CalcFEfromStrip(it->GetIndex()),
                                                     fStripCalcTop->CalcChannelfromStrip(it->GetIndex()),it->GetCharge());
        //dynamic_cast<PndHypDigiStrip*>(fStripArray->operator[](iStrip))->SetMCID(trackID);
        
        
        if (fVerbose > 1)std::cout << *it <<" "<<topStrips[iStrip] << std::endl;
        iStrip++;
      }
	  }else if(fVerbose>2) std::cout<<"Top side empty"<<std::endl;
	  
	  
	  // Bottom Side
	  //MVD x-y local coordinates,TVector2(0.,size.y()*2.)
	  //HYP x-z,y-z local coor. TVector2(0.,size.z()*2.)
	  if (fVerbose > 0)
	    std::cout  << "Bottom Side: " << std::endl;
    //  StripCalc = PndHypCalcStrip(fbotPitch, forient + fskew, fnrBotFE * fnrFECh , 
	  // 				   //FairGeoVector(-size.x(),-size.y(),0.),
	  // 				   TVector2(0.,size.y()*2.),
	  // 				   fthreshold, fnoise);
	  
	  fStripCalcBot->SetVerboseLevel(fVerbose);
	  Double_t stbot;
	  stbot = ((size.y()*2)/fDigiPar->GetBotPitch());
	  
	  Int_t rfebot = Int_t(stbot);Int_t rbfe = Int_t(rfebot/128);
    
	  if(rfebot>=rbfe*128)fStripCalcBot->SetreStrip(rbfe*128);
	  else fStripCalcBot->SetreStrip(rfebot);
    
	  std::cout<<" bof fe "<<Int_t(rfebot/128)<<" total strips "<<rfebot<<std::endl;
	  fStripCalcBot->SetAnchor(anchor);
	  
	  // Calculate a cluster of Strips fired
	  
	  std::vector<PndHypStrip> botStrips =
    fStripCalcBot->GetStrips(posInL.X(), posInL.Y(), posInL.Z(),
                             posOutL.X(), posOutL.Y(), posOutL.Z(),
                             point->GetEnergyLoss(),nSiL);
	  
	  if (botStrips.size() != 0){
	    //std::cout << "Deposited charge below threshold" << std::endl;
	    // continue;
	    //}
      // stripHits += botStrips.size();
      if (fVerbose > 1)
        std::cout  << "SensorStrips: bot " << botStrips.size()<<std::endl;
      
      trackID = point->GetTrackID();
      detID   = point->GetVolumeID();
      //uint iStrip = 0;
      Int_t sp;SensorSide si;
      
      for(std::vector<PndHypStrip>::const_iterator it=botStrips.begin();
          it!= botStrips.end();++it)
	    {
        new ((*fStripArray)[iStrip]) PndHypDigiStrip(iPoint,detID,
                                                     point->GetDetName(),
                                                     fStripCalcBot->CalcFEfromStrip(it->GetIndex()) + rtfe,
                                                     //StripCalcBot.CalcFEfromStrip(it->GetIndex()) + fnrTopFE,
                                                     fStripCalcBot->CalcChannelfromStrip(it->GetIndex()),it->GetCharge());
	      //    dynamic_cast<PndHypDigiStrip*>(fStripArray->operator[](iStrip))->SetMCID(trackID);
	      
        
        if (fVerbose > 1)
          std::cout << *it << std::endl;
	      
        iStrip++;
	    }
      
	  }else if(fVerbose>2) std::cout<<"Bottom side empty"<<std::endl;
	  // TODO read Parameters from Database
	} // Loop over MCPoints
  
  // Event summary
  std::cout << "-I- PndHypStripHitProducer: " << nPoints << " HypPoints, "
  <<" Hits created." << " " << iStrip<< std::endl;
  
}


// void PndHypStripHitProducer::GetLocalHitPoints(PndHypPoint* myPoint, FairGeoVector& myHitIn, FairGeoVector& myHitOut)
// {

//   if (fVerbose > 1)
//     std::cout << "GetLocalHitPoints" << std::endl;
//   //TGeoHMatrix trans = GetTransformation(myPoint->GetDetName().Data());
//   TGeoHMatrix trans = GetTransformation(fGeoH->GetPath(myPoint->GetDetName()).Data());

//   double posIn[3];
//   double posOut[3];
//   double posInLocal[3];
//   double posOutLocal[3];

//   posIn[0] = myPoint->GetXin();
//   posIn[1] = myPoint->GetYin();
//   posIn[2] = myPoint->GetZin();

//   posOut[0] = myPoint->GetXout();
//   posOut[1] = myPoint->GetYout();
//   posOut[2] = myPoint->GetZout();

//   if (fVerbose > 1){
//     for (int i = 0; i < 3; i++)
//       std::cout << "posIn "<< i << ": " << posIn[i] << std::endl;

//     trans.Print("");
//   }

//   trans.MasterToLocal(posIn, posInLocal);
//   trans.MasterToLocal(posOut, posOutLocal);

//   if (fVerbose > 1) {
//     for (int i = 0; i < 3; i++){
//       std::cout << "posInLocal "<< i << ": " << posInLocal[i] << std::endl;
//       std::cout << "posOutLocal "<< i << ": " << posOutLocal[i] << std::endl;
//     }
//   }


//   //posIn/OutLocal have the center of the coordinate system in the center of the shape
//   //typically sensors have their coordinate system centered at the lower left corner

//   TVector3 offset = GetSensorDimensions(fGeoH->GetPath(myPoint->GetDetName()).Data());
//   //TVector3 offset = GetSensorDimensions(myPoint->GetDetName().Data());

//   posInLocal[0] += offset.x();
//   posInLocal[1] += offset.y();
//   //posInLocal[2] += offset.z();

//   posOutLocal[0] += offset.x();
//   posOutLocal[1] += offset.y();
//   //posOutLocal[2] += offset.z();


//   myHitIn.setVector(posInLocal);
//   myHitOut.setVector(posOutLocal);

// }




// TGeoHMatrix PndHypStripHitProducer::GetTransformation(std::string detName) const
// {
//   gGeoManager->cd(detName.c_str());
//   TGeoHMatrix* transMat = gGeoManager->GetCurrentMatrix();
//   if (fVerbose > 1)
//     transMat->Print("");
//   return *transMat;
// }

TVector3 PndHypStripHitProducer::GetSensorDimensions(std::string detName) const
{
  gGeoManager->cd(detName.c_str());
  TGeoVolume* actVolume = gGeoManager->GetCurrentVolume();
  TGeoBBox* actBox = (TGeoBBox*)(actVolume->GetShape());
  TVector3 result;
  result.SetX(actBox->GetDX());
  result.SetY(actBox->GetDY());
  result.SetZ(actBox->GetDZ());
  
  //   //result.Dump();
  
  return result;
}
// -------------------------------------------------------------------------


ClassImp(PndHypStripHitProducer)