#include <stdio.h>
#include <cometEventLoop.hxx>
#include <TApplication.h>
#include <memory>
#include <TCanvas.h>
#include <TF1.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TBox.h>
#include <TGraph.h>
#include <TLine.h>
#include <TText.h>
#include <TMarker.h>
#include <TStyle.h>
#include "IMCHit.hxx"
#include <TFile.h>
#include <IOARuntimeParameters.hxx>
#include <IG4Trajectory.hxx>
#include <IG4PrimaryVertex.hxx>
#include <IG4TrajectoryPoint.hxx>

#include <IMultiHit.hxx>
#include <IValidationUtils.hxx>

#include <ICalibManager.hxx>

#include "IGlobalRecon.hxx"
#include <IFGDChannelCalibrator.hxx>
#include <IEventFolder.hxx>
#include <ICOMETInput.hxx>

#include <ReconPrintout.hxx>
#include <ReconObjectUtils.hxx>
#include <TrackingUtils.hxx>

TApplication *theApp;



TCanvas *WaveformCanvas;


#include <TGClient.h> 
#include <TGButton.h> 
#include <TGNumberEntry.h> 
#include <TGFrame.h> 

const int NSMRDVOLUMES=6;
const int NECALVOLUMES=13;

class IMyMainFrame : public TGMainFrame { 
 private: 
  TGNumberEntry *fNumber;
  TGTextButton  *fButton1;
  TGLayoutHints *fLayout;
  int val;
  
  int fRequestedEvent;

 public: 
  IMyMainFrame(const TGWindow *p, UInt_t w, UInt_t h); 
  ~IMyMainFrame(); 
  Bool_t ProcessMessage(Long_t msg, Long_t parm1, Long_t parm2); 
  void GoToEventNumber();

  int GetRequestedEvent(){return fRequestedEvent;}

};

IMyMainFrame::~IMyMainFrame(void)
{
  std::cout << " Entering delete Main Frame " << std::endl;
  delete fButton1;
  delete fNumber;

  delete fLayout;
}



//**************************************************************
class IBoxPad:public TBox {
public: 
  int tpc;
  int half;
  int mm; 
  int pad; 
  vector<double> time; 
  vector<double> adc;
  double timemin; 
  double timemax; 
  TH1F *histo; 

  IBoxPad();

  IBoxPad(COMET::IGeometryId geomId,Double_t x1, Double_t y1, Double_t x2, Double_t y2):TBox(x1,y1,x2,y2) {
    mm = COMET::IGeomInfo::Get().TPC().GeomIdToMM(geomId);
    pad  = COMET::IGeomInfo::Get().TPC().GeomIdToPad(geomId);
    tpc  = COMET::IGeomInfo::Get().TPC().GeomIdToTPC(geomId);
    half  = COMET::IGeomInfo::Get().TPC().GeomIdToHalf(geomId);
    timemin = 999999999;
    timemax = -timemin; 
    histo = NULL; 
  }


   
  ~IBoxPad(){
    if( histo ) 
      delete histo;

    histo = NULL;
 
  }

//**************************************************************
  void ExecuteEvent(Int_t event, Int_t px, Int_t py) {
//**************************************************************
    if (!gPad) return;
    
    TVirtualPad *savedpad = gPad; 

    if ( event != kMouseEnter) return;

    HideToolTip(event);
    
    TH1F *h =  histogram();
    if( h ) {
      WaveformCanvas->cd(); 
      h->Draw();
      WaveformCanvas->Update();
    }

    savedpad->cd();
  }


//**************************************************************
  void Add(COMET::IHandle<COMET::IHit> cHit){
//**************************************************************

    typedef std::vector< COMET::IHandle< COMET::ISingleHit > > Hits;
    typedef Hits::const_iterator iterator;
    COMET::IMultiHit *mh = dynamic_cast<COMET::IMultiHit *> (&(*cHit));
    
    if( mh ) {
      for( iterator mhit = mh->begin(); mhit !=  mh->end(); mhit++ ) {
  Add((*mhit)->GetTime(),(*mhit)->GetCharge());
      }
    }
    else {
      Add((cHit)->GetTime(),(cHit)->GetCharge());
    }
      
  }

//**************************************************************
  void Add(double t,double a) {
//**************************************************************
 
    if( t < timemin ) 
      timemin = t; 
    if( t > timemax ) 
      timemax = t;
    
    time.push_back(t);
    adc.push_back(a);
  }
  

//**************************************************************  
  TH1F *histogram(  ){
//**************************************************************

    if( timemin > timemax ) return NULL;

    if ( !histo ) {
      char name[256];
      char title[256];
      sprintf(name,"WF_%d_%d_%d_%d",tpc,half,mm,pad);
      sprintf(title,"TPC %d HALF %d MM %d PAD %d",tpc,half,mm,pad);

      histo = new TH1F(name,title,(int)(timemax-timemin+1),timemin,timemax);
      for(unsigned int i=0;i< time.size();i++){
  double a = adc[i];
  double t = time[i];
  histo->Fill(t,a);
      }
    }

    return histo; 
  }

};

 

vector<TBox*> MMvectorYZ;
vector<TBox*> MMvectorXZ;
vector<TBox*> MMvectorXY;
map<int,IBoxPad*> waveform; 

void DynamicExec()
{
  std::cout << " GOOD TRY " << std::endl; 
}



//**************************************************************
IMyMainFrame::IMyMainFrame(const TGWindow *p, UInt_t w, UInt_t h):TGMainFrame(p,w,h){
//**************************************************************

// Create a main frame with a number of different buttons.
  fButton1 = new TGTextButton(this, "&Next", 1);
  fNumber = new TGNumberEntry(this, 999, 7,999, TGNumberFormat::kNESInteger,
            TGNumberFormat::kNEANonNegative,TGNumberFormat::kNELLimitMinMax,0, 999999);

  fLayout = new TGLayoutHints(kLHintsCenterX | kLHintsCenterY | kLHintsExpandX | kLHintsExpandY );
  AddFrame(fButton1, fLayout);
  AddFrame(fNumber, fLayout);
  MapSubwindows();
  Layout();
  SetWindowName("Button Control");
  SetIconName("Button Control ");
  MapWindow();
}


//**************************************************************
Bool_t IMyMainFrame::ProcessMessage(Long_t msg, Long_t parm1, Long_t){
//**************************************************************

  // Process events generated by the buttons in the frame.
  switch (GET_MSG(msg)){
  case kC_COMMAND:
    switch (GET_SUBMSG(msg)){
    case kCM_BUTTON:
      printf("text button id %ld pressed\n", parm1);
      if( parm1 == 1 ) {
  fRequestedEvent=0;
  std::cout << " New Event " << std::endl; 
  theApp->Terminate(true);
      }      
    default:
      break;
    }
  case kC_TEXTENTRY:
    switch (GET_SUBMSG(msg)) {
    case kTE_ENTER:
      fRequestedEvent=fNumber->GetIntNumber();
      std::cout << "requestedEvent: " << fRequestedEvent << std::endl;
      theApp->Terminate(true);
    }	  
  default:
    break;

  }

  gROOT->Reset();

 return kTRUE;
}

//**************************************************************
void IMyMainFrame::GoToEventNumber(){
//**************************************************************

  Int_t requestedEvent=fNumber->GetIntNumber();
  std::cout << "requestedEvent: " << requestedEvent << std::endl;

}


//**************************************************************
class IRunReconGlobalEvDisplay: public COMET::ICOMETEventLoopFunction{
//**************************************************************

public:
  IRunReconGlobalEvDisplay();
  IRunReconGlobalEvDisplay(const TGWindow *p, UInt_t w, UInt_t h);
  
  virtual ~IRunReconGlobalEvDisplay() {}

  void GeomInit();

  bool isMC(COMET::ICOMETEvent& event);

  void DrawDetector();
  void DrawEventInfo( COMET::ICOMETEvent& event);
  void DrawTPCHits(COMET::IHandle<COMET::IHitSelection> tpc, double t0, int hitColor);
  void Draw2DHits(COMET::IHandle<COMET::IHitSelection> hits, int hitColor);
  void DrawNodesInObject(COMET::IHandle<COMET::IReconBase> object, int nodeColor);
  void DrawMCTrajectories(COMET::IHandle<COMET::IReconBase> object, int truthColor);
  void DrawMCVertices(COMET::ICOMETEvent& event, int vertexColor);
  void DrawRecVertices(COMET::IHandle<COMET::IReconObjectContainer> vertices, int vertexColor);
  void DrawVertices(int npoints, double *x, double *y, double *z, int color, int marker);
  void DrawPoints(int npoints, double *x, double *y, double *z, int color);

  bool FindEvent(int requestedEvent, COMET::ICOMETEvent* event0, COMET::ICOMETEvent*& event);

  void InitExplanation();
  void DrawExplanation();

  COMET::IHandle<COMET::IHitSelection> GetTPCHits(COMET::ICOMETEvent& event);
  COMET::IHandle<COMET::IHitSelection> GetFGDHits(COMET::ICOMETEvent& event);
  COMET::IHandle<COMET::IHitSelection> GetHitsFromComboHits(COMET::IHandle<COMET::IHitSelection> tpc);

  void BeginFile(COMET::IVInputFile *input);
  void Initialized();
  bool SetOption(std::string option, std::string value);

  bool operator () (COMET::ICOMETEvent& event);

  IMyMainFrame *mf;
  
private:
  int NBINS;

  IGlobalRecon* fGlobalRecon;

  TH2F* YZHisto;
  TH2F* XZHisto;
  TH2F* XYHisto;

  TCanvas *YZ;
  TCanvas *XZ;
  TCanvas *XY;


  std::vector<TBox*> box_XZ;
  std::vector<TBox*> box_YZ;
  std::vector<TBox*> box_XY;


  TLine *tline[5];
  TText *text[7];
  TGraph *tmarker[2];

  TText *EvtLabel;

  vector<IBoxPad*> fBoxPadList;
  vector<TGraph*>  fGraphList;

  bool fForceRunGlobalRecon; 
  bool fForceRunSubdetectorRecon; 

  double fCathodeOffset; 

  double fTimeOffset;

  double fDriftVelocity;
  double fMaxDrift;
  double fSampleThreshold;

  double fSamplingTime;
  

  int fUnusedHitsColor;
  int fOtherHitsColor;
  int fObjectHitsColor;
  int fObjectNodesColor;
  int fTruthColor;
  int fTrueVertexColor;
  int fRecVertexColor;

  std::string fRecPackVerb;
  int ftmanVerb;
  int fConverterVerb;
  int fpmanVerb;

  bool geomInitialized;

  COMET::ICOMETInput* fInputFile;
  
};


//*************************************************
IRunReconGlobalEvDisplay::IRunReconGlobalEvDisplay() {
//*************************************************

  COMET::ICOMETLog::SetLogLevel(COMET::ICOMETLog::QuietLevel);
  COMET::ICOMETLog::SetLogLevel("ReconGlobal",COMET::ICOMETLog::QuietLevel);    

  fRecPackVerb ="00000";
  ftmanVerb = 0;
  fConverterVerb = 0;
  fpmanVerb = 0;

  COMET::ICalibManager::Get().AddCalibrator(new COMET::IFGDChannelCalibrator());

  geomInitialized = false;

}

//*************************************************
bool IRunReconGlobalEvDisplay::SetOption(std::string option, std::string value) {
//*************************************************

  COMET::ICOMETLog::LogPriority verb[4] = {COMET::ICOMETLog::QuietLevel, COMET::ICOMETLog::LogLevel, COMET::ICOMETLog::InfoLevel, COMET::ICOMETLog::VerboseLevel};

  int value2 = atoi(value.c_str());
  if (value2 <1) value2=0;
  if (value2 >3) value2=3;

  if (option == "v"){
    COMET::ICOMETLog::SetLogLevel("ReconGlobal",verb[value2]);
  }
  else if (option == "vorp"){
    COMET::ICOMETLog::SetLogLevel("oaRecPack",verb[value2]);
  }
  else if (option == "vtr"){
    COMET::ICOMETLog::SetLogLevel("ReconTracker",verb[value2]);
  }
  else if (option == "vfgd"){
    COMET::ICOMETLog::SetLogLevel("ReconFGD",verb[value2]);
  }
  else if (option == "vsmrd"){
    COMET::ICOMETLog::SetLogLevel("ReconSMRD",verb[value2]);
  }
  else if (option == "vrp"){
    fRecPackVerb = value;
  }
  else if (option == "vtman"){
    ftmanVerb = atoi(value.c_str());
  }
  else if (option == "vconv"){
    fConverterVerb = atoi(value.c_str());
  }
  else if (option == "vpman"){
    fpmanVerb = atoi(value.c_str());
  }    
  else return false;

  return true;

}


//*****************************************************************************
void IRunReconGlobalEvDisplay::BeginFile (COMET::IVInputFile *input){
//*****************************************************************************

  fInputFile = dynamic_cast<COMET::ICOMETInput*>(input);
    
}

//*****************************************************************************
void IRunReconGlobalEvDisplay::GeomInit (){
//*****************************************************************************


  fGlobalRecon = new IGlobalRecon();


  // dums the RecPack Setup
  if(COMET::ICOMETLog::GetLogLevel("ReconGlobal") >= COMET::ICOMETLog::VerboseLevel){
    std::cout << COMET::tman().GetSetup() << std::endl;  
    COMET::gman().DumpVolumeProperties();
  }

  // Changes the RecPack verbosity (5 digits = fitting, navigation, model, matching, RayTool. Minimum is 0, maximum is 7)
  COMET::tman("globalRecon").SetVerbosity(fRecPackVerb);   

  if(COMET::ICOMETLog::GetLogLevel("ReconTracker") >= COMET::ICOMETLog::VerboseLevel)
    COMET::tman("ReconTracker").SetVerbosity(fRecPackVerb);   

  if(COMET::ICOMETLog::GetLogLevel("ReconFGD") >= COMET::ICOMETLog::VerboseLevel)
    COMET::tman("ReconFGD").SetVerbosity(fRecPackVerb);   

  if(COMET::ICOMETLog::GetLogLevel("ReconSMRD") >= COMET::ICOMETLog::VerboseLevel)
    COMET::tman("ReconSMRD").SetVerbosity(fRecPackVerb);   
 
  // the IRecPackConverters verbosity
  COMET::converter().SetVerbosity(fConverterVerb);    
  
  // the IPIDManager verbosity
  COMET::pman().SetVerbosity(fpmanVerb);    


  const TGeoNode *cath = COMET::IGeomInfo::Get().TPC().GetCathode(0);
  const TGeoVolume *cathvol =	cath->GetVolume();
  TGeoShape       *cathshape =	cathvol->GetShape();
  TGeoBBox* box = dynamic_cast<TGeoBBox*>(cathshape);
  fCathodeOffset =  box->GetDX(); 
  
  std::cout << " CATHODE " << fCathodeOffset << std::endl; 


  WaveformCanvas =  new TCanvas ("WaveForm","  WaveForm  ",520,420,400,400);

  
  YZ = new TCanvas ("YZ","  YZ  ",100,0,600,400);
  XZ = new TCanvas ("XZ","  XZ  ",720,0,600,400);
  XY = new TCanvas ("XY","  XY  ",100,420,400,400);
  
  YZHisto = new TH2F("YZHisto"," YZ Projection ",100,-4000.,4000.,100,-2400.,3400.);
  XZHisto = new TH2F("XZHisto"," XZ Projection ",100,-4000.,4000.,100,-2400.,2400.);
  XYHisto = new TH2F("XYHisto"," XY Projection ",100,-2400.,2400.,100,-2600.,2600.);  

}



//**************************************************************  
void IRunReconGlobalEvDisplay::Initialized() {
//**************************************************************

  int argc2 = 0;
  char **argv2 = NULL;

  fGlobalRecon = 0;
  
  theApp = new TApplication("App",&argc2,argv2);
  
  gROOT->SetStyle("Plain");
  gStyle->SetOptStat(0);
  
  mf = new IMyMainFrame(gClient->GetRoot(),100, 100);
  
  mf->SetWMPosition(10,10);

  
  fCathodeOffset = 0.0;

  fForceRunGlobalRecon      = (bool)COMET::IOARuntimeParameters::Get().GetParameterI("ReconGlobal.ForceRunGlobalRecon");
  fForceRunSubdetectorRecon = (bool)COMET::IOARuntimeParameters::Get().GetParameterI("ReconGlobal.ForceRunSubdetectorRecon");

  fTimeOffset = COMET::IOARuntimeParameters::Get().GetParameterD("ReconTPC.Reco.LikFit.TimeOffset");
  fDriftVelocity =  COMET::IOARuntimeParameters::Get().GetParameterD("ReconTPC.TPCgas.DriftSpeed");
  fMaxDrift     =   COMET::IGeomInfo::Get().TPC().GetMaxDriftDistance(0,0);
  fSampleThreshold = COMET::IOARuntimeParameters::Get().GetParameterD("ReconTPC.Reco.Cluster.SampleThreshold");

  fSamplingTime = COMET::IOARuntimeParameters::Get().GetParameterD("ReconTPC.AfterTPC.SamplingTime");


  fUnusedHitsColor = 2;
  fOtherHitsColor = 6;
  fObjectHitsColor = 3;
  fObjectNodesColor = 4;
  //  fTruthColor = 5;
  fTruthColor = 41;
  fTrueVertexColor = 41;
  fRecVertexColor = 1;


  EvtLabel = NULL; 

  InitExplanation();


  std::cout << "fTimeOffset:    "<< fTimeOffset << std::endl;
  std::cout << "fDriftVelocity: "<< fDriftVelocity << std::endl;
  std::cout << "fMaxDrift:      "<< fMaxDrift << std::endl;
  
}


//**************************************************************  
void IRunReconGlobalEvDisplay::DrawDetector() {
//**************************************************************  


  // -------------- Build the TPC MM ------------------    
  std::map<int,std::pair<int,const TGeoNode*> >  MM = COMET::IGeomInfo::Get().TPC().GetMMVolumes();
  
  if( MMvectorYZ.size() == 0 ) {
    
    std::cout << " MM modules " << MM.size() << std::endl; 
    
    int saved = gGeoManager->GetNodeId();
    
    for(std::map<int,std::pair<int,const TGeoNode*> >::iterator ib = MM.begin(); ib != MM.end();ib++){
      
      const TGeoNode *mm = (ib->second).second;
      const TGeoVolume *mvol =	mm->GetVolume();
      TGeoShape       *shape =	mvol->GetShape();
      TGeoBBox* box = dynamic_cast<TGeoBBox*>(shape);
      double x = box->GetDX();
      double y = box->GetDY();
      double z = box->GetDZ();
      Double_t local[3],master[3];
      
      local[0] = local[1] = local[2] =  0.0;
      
      gGeoManager->CdNode((ib->second).first);
      
      gGeoManager->LocalToMaster(local,master);
      
      TBox *nn;
      if (master[0] < 0){       
  nn = new TBox(master[2]-x,master[1]-y,master[2]+x,master[1]+y);
  
  MMvectorYZ.push_back(nn);
  nn->SetFillStyle(0);
  nn->SetLineColor(1);
      }	  
      
      nn = new TBox(master[2]-x,master[0]-z,master[2]+x,master[0]+z);
      
      MMvectorXZ.push_back(nn);
      nn->SetFillStyle(0);
      nn->SetLineColor(1);
      
      
      nn = new TBox(master[0]-z,master[1]-y,master[0]+z,master[1]+y);
      
      MMvectorXY.push_back(nn);
      nn->SetFillStyle(0);
      nn->SetLineColor(1);
      
    }
    
    gGeoManager->CdNode(saved);    

  }
  
  const std::vector<dict::Key>&  volumes = COMET::gman().GetSetup().volumes().keys(); 
  box_XZ.resize(volumes.size());
  box_XY.resize(volumes.size());
  box_YZ.resize(volumes.size());

  for (unsigned int i=0;i<volumes.size();i++){

    if (volumes[i]=="mother" || volumes[i]=="/t2k_1/OA_0/Magnet_0") continue;

    const Volume& vol = COMET::gman().GetSetup().volume(volumes[i]);
    const EVector size = vol.size();
    const EVector pos = vol.position();
    TVector3 smin(pos[0]-size[0],pos[1]-size[1],pos[2]-size[2]);
    TVector3 smax(pos[0]+size[0],pos[1]+size[1],pos[2]+size[2]);
    box_XZ[i] = new TBox(smin.Z(),smin.X(), smax.Z(), smax.X());
    box_YZ[i] = new TBox(smin.Z(),smin.Y(), smax.Z(), smax.Y());
    box_XY[i] = new TBox(smin.X(),smin.Y(), smax.X(), smax.Y());    

    int color=1;
    int fcolor=0;
    int fcolor2=0;
    int fstyle=0;
    int fstyle2=0;

    if (volumes[i].find("TPC") !=std::string::npos ) color=1;
    if (volumes[i].find("FGD") !=std::string::npos ) color=2;
    if (volumes[i].find("DsECal") !=std::string::npos ) color=3;
    if (volumes[i].find("P0DECal") !=std::string::npos ) color=4;
    if (volumes[i].find("BrlECal") !=std::string::npos ) color=8;
    if (volumes[i].find("P0D_0") !=std::string::npos ) color=6;
    if (volumes[i].find("SMRD") !=std::string::npos ) color=7;

    fstyle=0;
    fcolor=color;

    fstyle2=fstyle;
    fcolor2=fcolor;

    /*
    if (volumes[i].find("Active") !=std::string::npos ||
  ( volumes[i].find("Drift") !=std::string::npos && volumes[i].find("Cathode") ==std::string::npos)  ){
  fcolor2=18;
  fstyle2=1001;
      }
    */
    box_XZ[i]->SetFillStyle(fstyle2);
    box_YZ[i]->SetFillStyle(fstyle2);
    box_XY[i]->SetFillStyle(fstyle);

    box_XZ[i]->SetFillColor(fcolor2);
    box_YZ[i]->SetFillColor(fcolor2);
    box_XY[i]->SetFillColor(fcolor);

    box_XZ[i]->SetLineColor(color);
    box_YZ[i]->SetLineColor(color);
    box_XY[i]->SetLineColor(color);

    box_XZ[i]->SetLineStyle(1);
    box_YZ[i]->SetLineStyle(1);
    box_XY[i]->SetLineStyle(1);


  }


  YZHisto->GetXaxis()->UnZoom();
  YZHisto->GetYaxis()->UnZoom();
  XYHisto->GetYaxis()->UnZoom();
  XYHisto->GetXaxis()->UnZoom();
  XZHisto->GetXaxis()->UnZoom();
  XZHisto->GetYaxis()->UnZoom();


  YZ->cd();
  YZHisto->Draw();
  for (unsigned int i=0;i<volumes.size();i++){
    if (volumes[i].find("Side") != std::string::npos || 
  (volumes[i].find("Left") != std::string::npos && volumes[i].find("LeftClam")  == std::string::npos)|| 
  (volumes[i].find("Right")!= std::string::npos && volumes[i].find("RightClam") == std::string::npos)) continue;
    if (box_YZ[i]) box_YZ[i]->Draw("same");
  }

  for(unsigned int i = 0;i< MMvectorYZ.size() ; i++) 
    MMvectorYZ[i]->Draw("same");


  YZ->Update();

  XZ->cd();
  XZHisto->Draw();
  for (unsigned int i=0;i<volumes.size();i++){
    if (volumes[i].find("Top") !=std::string::npos ) continue;
    if (volumes[i].find("Bottom") !=std::string::npos ) continue;
    if (box_XZ[i]) box_XZ[i]->Draw("same");
  }

  for(unsigned int i = 0;i< MMvectorXZ.size() ; i++) 
    MMvectorXZ[i]->Draw("same");


  XZ->Update();

  XY->cd();
  XYHisto->Draw();
  for (unsigned int i=0;i<volumes.size();i++)
    if (box_XY[i]) box_XY[i]->Draw("same");

  for(unsigned int i = 0;i< MMvectorXY.size() ; i++) 
    MMvectorXY[i]->Draw("same");


  XY->Update();

  return;


  
}



//**************************************************************
bool IRunReconGlobalEvDisplay::isMC(COMET::ICOMETEvent& event) {
//**************************************************************

  //    COMET::IHandle<COMET::IG4TrajectoryContainer> g4traj = event.Get<COMET::IG4TrajectoryContainer>("truth/G4Trajectories");
  COMET::IHandle<COMET::IG4PrimaryVertexContainer> g4vert
    = event.Get<COMET::IG4PrimaryVertexContainer>("truth/G4PrimVertex00");
  
  if( g4vert) return true;
  
  return false;
}


//**************************************************************
void IRunReconGlobalEvDisplay::DrawEventInfo( COMET::ICOMETEvent& event){
//**************************************************************
  
  int trigbits =(int) event.GetHeader().GetTriggerBits();
  int run = event.GetRunId();
  int subrun = event.GetContext().GetSubRun();
  int evt = event.GetEventId();
  
  //  fRun = run; 
  //  fEvent = evt;
  
  char label[256]; 
  
  if (subrun<0) subrun=0;

  sprintf(label," Run %d, SRun %d, Evt %d, Trig 0x%x",run,subrun, evt,trigbits);
  
  if( EvtLabel ) delete EvtLabel; 
  
  EvtLabel = new TText(-3800.,2600.,label);
  
  YZ->cd();
  EvtLabel->Draw(); 
  YZ->Update();

}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawTPCHits(COMET::IHandle<COMET::IHitSelection> tpc, double t0, int hitColor) {
//**************************************************************

  if (!tpc) return;

  COMET::IHandle<COMET::IHitSelection> hits =  GetHitsFromComboHits(tpc);  

  map<int,bool> used;   
  for(COMET::IHitSelection::const_iterator hsel = hits->begin(); hsel != hits->end(); hsel++){

    if (TrackingUtils::GetDetector(*hsel) != COMET::IReconBase::kTPC1 &&
  TrackingUtils::GetDetector(*hsel) != COMET::IReconBase::kTPC2 &&
  TrackingUtils::GetDetector(*hsel) != COMET::IReconBase::kTPC3) continue;

    COMET::IGeometryId geomId = (*hsel)->GetGeomId();
    const TVector3 pos = (*hsel)->GetPosition();
    IBoxPad *tb;
    int localid = geomId.AsInt();
    
    if( (*hsel)->GetCharge() < fSampleThreshold ) continue;  // For data mainly, reduce the number of plotted hits to a handlable amount. 

    // YZ projection
    if( !used[localid] ){
      used[localid] = true;
      YZ->cd();
      tb = new IBoxPad( geomId,pos.Z()-9.7/2., pos.Y()-6.9/2.,pos.Z()+9.7/2., pos.Y()+6.9/2.);
      fBoxPadList.push_back(tb);
      tb->SetFillStyle(1);
      tb->SetLineColor(hitColor);
      tb->SetFillColor(hitColor);
      tb->Draw("same");
      waveform[localid] = tb;
      waveform[localid]->Add(*hsel);
    } else
      waveform[localid]->Add(*hsel);


    double peakTime=0;
    double maxCharge = 0; 
    
    // Find the pick charge and Time of the Wave form
    COMET::IHandle<COMET::IMultiHit> mh = *hsel;
    std::vector< COMET::IHandle< COMET::ISingleHit > >::const_iterator mhit;
    for( mhit = mh->begin(); mhit !=  mh->end(); mhit++ ) {
      if( maxCharge < (*mhit)->GetCharge() ) {
  maxCharge = (*mhit)->GetCharge(); 
  peakTime = (*mhit)->GetTime();
      }
    }
    
    double t = peakTime - t0;
    double x = (fMaxDrift + fCathodeOffset - t*fDriftVelocity)*COMET::IGeomInfo::Get().TPC().GetDriftSense(geomId);


    // XZ projection
    XZ->cd();
    TVector3 global;
    COMET::IGeomInfo::Get().TPC().GeomIdToGlobalXYZ(geomId,global);
    tb = new IBoxPad( geomId, pos.Z()-9.7/2.,x - 50.*fDriftVelocity/2.,pos.Z()+9.7/2.,x + 50.*fDriftVelocity/2.);
    fBoxPadList.push_back(tb);
    tb->SetFillStyle(1);
    tb->SetLineColor(1);
    tb->SetFillColor(hitColor);
    tb->Draw("same");
    
    // XY projection
    XY->cd();
    COMET::IGeomInfo::Get().TPC().GeomIdToGlobalXYZ(geomId,global);
    tb = new IBoxPad( geomId, x - 50.*fDriftVelocity/2.,pos.Y()-6.8/2.,x + 50.*fDriftVelocity/2.,pos.Y()+6.8/2.);
    fBoxPadList.push_back(tb);
    tb->SetFillStyle(1);
    tb->SetLineColor(hitColor);
    tb->SetFillColor(hitColor);
    tb->Draw("same");
    
  }
}


//**************************************************************
void IRunReconGlobalEvDisplay::Draw2DHits(COMET::IHandle<COMET::IHitSelection> hits, int hitColor) {
//**************************************************************

  
  //------------ Draw reconstructed 2D hits ------------------------------      
  if( !hits ) return;
    
  for(COMET::IHitSelection::const_iterator hsel = hits->begin(); hsel != hits->end(); hsel++){

    if (TrackingUtils::GetDetector(*hsel) == COMET::IReconBase::kTPC1 ||
  TrackingUtils::GetDetector(*hsel) == COMET::IReconBase::kTPC2 ||
  TrackingUtils::GetDetector(*hsel) == COMET::IReconBase::kTPC3) continue;

    const TVector3 pos = (*hsel)->GetPosition();

    if  ((*hsel)->GetCharge()<2) continue;
    
    if( (*hsel)->IsYHit() ){
      YZ->cd();
      IBoxPad *tb = new IBoxPad( COMET::IGeometryId(0), pos.Z()-10.0/2., pos.Y()-10.0/2., pos.Z()+10.0/2., pos.Y()+10.0/2.);
      fBoxPadList.push_back(tb);
      tb->SetFillStyle(1);
      tb->SetLineColor(1);
      tb->SetFillColor(hitColor);
      tb->Draw("same");
    }
    else if( (*hsel)->IsXHit()     ){
      
      XZ->cd();	    
      IBoxPad *tb = new IBoxPad( COMET::IGeometryId(0),pos.Z()-10.0/2., pos.X()-10.0/2.,pos.Z()+10.0/2., pos.X()+10.0/2.);
      fBoxPadList.push_back(tb);
      tb->SetFillStyle(1);
      tb->SetLineColor(1);
      tb->SetFillColor(hitColor);
      tb->Draw("same");
    }
  }
  
}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawNodesInObject(COMET::IHandle<COMET::IReconBase> object, int nodeColor){
//**************************************************************


  if(object->GetNodes().size()==0) return;


  Trajectory traj;
  COMET::converter().TReconBase_to_Trajectory(object, traj);
  //  std::cout << traj << std::endl;

  if (!traj.status("fitted")) return;

  COMET::IHandle<COMET::IReconState> newstate;

  COMET::IHandle<COMET::IReconPID> pid = object; 
  if (pid){
    //    COMET::converter().Trajectory_to_TReconBase(traj, "PID", object);
     newstate= COMET::IHandle<COMET::IReconState>(new COMET::IPIDState);
  }
  else{
    COMET::IHandle<COMET::IReconTrack> track = object; 
    if (track){
      //      COMET::converter().Trajectory_to_TReconBase(traj, "Track", object);
      newstate= COMET::IHandle<COMET::IReconState>(new COMET::ITrackState);
    }
    else{ 
      COMET::IHandle<COMET::IReconShower> shower = object; 
      if (shower){
  //	COMET::converter().Trajectory_to_TReconBase(traj, "Shower", object);
  newstate= COMET::IHandle<COMET::IReconState>(new COMET::IShowerState);
      }
      else 		
  return;
    }
  }


  COMET::IReconNodeContainer &nodes = object->GetNodes();
  int npoints = nodes.size();

  double x[500];
  double y[500];
  double z[500];
  
  int i = 0;
  
  // Draw the nodes. 

  TVector3 normal(0,0,1);
  bool first_node=true;
  COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_length_sign(0);
  COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_allow_zero_length(true);
  for( COMET::IReconNodeContainer::iterator  ni = nodes.begin(); ni != nodes.end();ni++) {
    if (!((*ni)->GetState()) ) continue;

    COMET::IHandle<COMET::IReconBase> nodeObject = (*ni)->GetObject();

    //    std::cout << i << " node: " << nodeObject << std::endl;
    //    EventUtils::DumpReconHistory(nodeObject);
    //    if (nodeObject->GetAlgorithmName().find(RP::pos_slope_curv) != std::string::npos && 
    if (nodeObject->CheckStatus(COMET::IReconBase::kLikelihoodFit) && 
  nodeObject->UsesDetector(COMET::IReconBase::kTPC)){
      //      std::cout << i << ": " << COMET::tman().GetPosition((*ni)->GetState()).Vect().Z() << " " << COMET::tman().GetDetector(COMET::tman().GetPosition((*ni)->GetState()).Vect()) << " " << nodeObject->ConvertDetector() << std::endl;
      //      COMET::IHandle<COMET::IReconBase> object2 = COMET::tman().GetConstituentInDetector(object, COMET::tman().GetDetector(COMET::tman().GetPosition((*ni)->GetState()).Vect() ) );
      COMET::IHandle<COMET::IReconBase> object2 = *(nodeObject->GetConstituents()->begin());
      if (object2){
  //	std::cout << i << " node': " << object2 << std::endl;
  
  COMET::IHandle<COMET::IReconState> state = (*ni)->GetState();

  
  COMET::IReconNodeContainer &nodes2 = object2->GetNodes();
  if (first_node || !TrackingUtils::GetCurvBackTrack(*object)){
    first_node=false;
    //	  std::cout << " first node"<< std::endl;
    bool first_node2=true;
    for( COMET::IReconNodeContainer::iterator  ni2 = nodes2.begin(); ni2!= nodes2.end();ni2++) {
      if (!((*ni2)->GetState()) ) continue;
      
      TLorentzVector OldNodePos = TrackingUtils::GetPosition((*ni2)->GetState());  
      
      //	    std::cout << "old node pos: " << OldNodePos.X() << " "  << OldNodePos.Y() << " "  << OldNodePos.Z() << std::endl;
      
      TVector3 dir = TrackingUtils::GetDirection(state);
      if (first_node2){
        if (TrackingUtils::GetCurvBackTrack(*object))
    COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_length_sign(1);
        normal = TVector3(0,0,1);
        first_node2=false;
      }
      else{ 
        COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_length_sign(0);	  
        if (fabs(dir.Z()) > 0.3)
    normal = TVector3(0,0,1);
        else if (fabs(dir.Z()) > fabs(dir.Y()) &&  fabs(dir.Z()) > fabs(dir.X()))
    normal = TVector3(0,0,1);
        else if (fabs(dir.X()) > fabs(dir.Y()) &&  fabs(dir.X()) > fabs(dir.Z()))
    normal = TVector3(1,0,0);
        else if (fabs(dir.Y()) > fabs(dir.X()) &&  fabs(dir.Y()) > fabs(dir.Z()))
    normal = TVector3(0,1,0);
      }	    

      if(!COMET::tman().PropagateToSurface(*state,OldNodePos.Vect(),normal,*newstate)) continue;
      TLorentzVector nodePos = TrackingUtils::GetPosition(newstate);
      
      //	    std::cout << "new node pos: " << nodePos.X() << " "  << nodePos.Y() << " "  << nodePos.Z() << std::endl;
      
      x[i] = nodePos.X();
      y[i] = nodePos.Y();
      z[i] = nodePos.Z();
      
      i++;
      
      *state = *newstate;
    }
  }
  else{
    //	  std::cout << " second node"<< std::endl;
    bool first_node2=true;
    for( COMET::IReconNodeContainer::reverse_iterator  ni2 = nodes2.rbegin(); ni2!= nodes2.rend();ni2++) {
      if (!((*ni2)->GetState()) ) continue;
      
      TLorentzVector OldNodePos = TrackingUtils::GetPosition((*ni2)->GetState());  
      
      //	    std::cout << "old node pos: " << OldNodePos.X() << " "  << OldNodePos.Y() << " "  << OldNodePos.Z() << std::endl;
      
      TVector3 dir = TrackingUtils::GetDirection(state);
      if (first_node2){
        if (TrackingUtils::GetCurvBackTrack(*object))
    COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_length_sign(-1);
        normal = TVector3(0,1,0);
        first_node2=false;
      }
      else{
        COMET::rpman().model_svc().model(COMET::rpman().model_svc().model_name()).intersector().set_length_sign(0);
        if (fabs(dir.Z()) > 0.6)
    normal = TVector3(0,0,1);
        else if (fabs(dir.Z()) > fabs(dir.Y()) &&  fabs(dir.Z()) > fabs(dir.X()))
    normal = TVector3(0,0,1);
        else if (fabs(dir.X()) > fabs(dir.Y()) &&  fabs(dir.X()) > fabs(dir.Z()))
    normal = TVector3(1,0,0);
        else if (fabs(dir.Y()) > fabs(dir.X()) &&  fabs(dir.Y()) > fabs(dir.Z()))
    normal = TVector3(0,1,0);
      }	      

      if(!COMET::tman().PropagateToSurface(*state,OldNodePos.Vect(),normal,*newstate)) continue;
      TLorentzVector nodePos = TrackingUtils::GetPosition(newstate);
      
      //	    std::cout << "new node pos: " << nodePos.X() << " "  << nodePos.Y() << " "  << nodePos.Z() << std::endl;
      
      x[i] = nodePos.X();
      y[i] = nodePos.Y();
      z[i] = nodePos.Z();
      
      i++;
      
      *state = *newstate;
    }
  }
      }
      
    }
    else{
      first_node = true;

      TLorentzVector pos = TrackingUtils::GetPosition((*ni)->GetState());  
      x[i] = pos.X();
      y[i] = pos.Y();
      z[i] = pos.Z();
      
      i++;
    }
  }

  npoints = i;

  if (npoints==0) return;

  double xf[npoints];
  double yf[npoints];
  double zf[npoints];

  for (int j=0;j<npoints;j++){
    xf[j]=x[j];
    yf[j]=y[j];
    zf[j]=z[j];
    //    std::cout << j << ": " << xf[j] << " " << yf[j] << " " << zf[j] << std::endl;
  }

  DrawPoints(npoints,xf,yf,zf,nodeColor);
  
}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawMCTrajectories(COMET::IHandle<COMET::IReconBase> object, int truthColor) {
//**************************************************************

  
  double eff, pur;
  COMET::IHandle<COMET::IG4Trajectory> g4traj = COMET::vman().GetG4Trajectory(*object,pur,eff);
  if (!g4traj) return;

  double p = g4traj->GetInitialMomentum().Vect().Mag();
  TVector3 u = (1/p)*g4traj->GetInitialMomentum().Vect();
  TVector3 r = g4traj->GetInitialPosition().Vect();

  std::cout << object << std::endl;
  std::cout << "   --> Rec-Truth matching purity = " << pur << std::endl;
  std::cout << "      - True position =  " << r.X() << " " << r.Y() << " " << r.Z() << std::endl;
  std::cout << "      - True direction = " << u.X() << " " << u.Y() << " " << u.Z() << std::endl;
  std::cout << "      - True momentum =  " << p << std::endl;


  int npoints = g4traj->GetTrajectoryPoints().size();

  if (npoints==0) return;

  // Draw the MC true information if available.
  double x[npoints];
  double y[npoints];
  double z[npoints];
  
  std::vector<COMET::IG4TrajectoryPoint> pmc = g4traj->GetTrajectoryPoints();
  
  int i = 0;
  for( std::vector<COMET::IG4TrajectoryPoint>::iterator tp = pmc.begin();tp != pmc.end();tp++){
    TLorentzVector pos = tp->GetPosition();  
    x[i] = pos.X();
    y[i] = pos.Y();
    z[i] = pos.Z();
    
    i++;
  }

  DrawPoints(npoints,x,y,z,truthColor);
  
}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawMCVertices(COMET::ICOMETEvent& event, int vertexColor) {
//**************************************************************

  if( ! (event.Has<COMET::IG4PrimaryVertexContainer>("truth/G4PrimVertex00"))) return;
    
  std::cerr<<"Filling MC primary vertex plot..."<<std::endl;
  
  COMET::IHandle<COMET::IG4PrimaryVertexContainer>  g4PrimVert  = event.Get<COMET::IG4PrimaryVertexContainer>("truth/G4PrimVertex00");     
  int npoints = g4PrimVert->size();

  if (npoints==0) return;
  
  double x[npoints];
  double y[npoints];
  double z[npoints];
  
  for(unsigned int i = 0; i < (*g4PrimVert).size(); i++){
    TLorentzVector pos = (*g4PrimVert)[i].GetPosition();      
    x[i] = pos.X();
    y[i] = pos.Y();
    z[i] = pos.Z();
  }

  DrawVertices(npoints,x,y,z,vertexColor,kFullStar);
  
}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawRecVertices(COMET::IHandle<COMET::IReconObjectContainer> vertices, int vertexColor) {
//**************************************************************

  if (!vertices) return;

  std::cerr<<"Filling Rec vertex plot..."<<std::endl;
  
  int npoints = vertices->size();
  
  if (npoints==0) return;

  double x[npoints];
  double y[npoints];
  double z[npoints];
  
  for(unsigned int i = 0; i < vertices->size(); i++){
    COMET::IHandle<COMET::IReconVertex> vertex = (*vertices)[i]; 
    if (!vertex) continue;
    TLorentzVector pos = vertex->GetPosition();      
    x[i] = pos.X();
    y[i] = pos.Y();
    z[i] = pos.Z();
  }

  DrawVertices(npoints,x,y,z,vertexColor,kOpenCircle);

}


//**************************************************************
void IRunReconGlobalEvDisplay::DrawPoints(int npoints, double *x, double *y, double *z, int color) {
//**************************************************************
  
  if (npoints ==0) return;

  TGraph *txy = new TGraph(npoints,x,y);
  TGraph *txz = new TGraph(npoints,z,x);
  TGraph *tyz = new TGraph(npoints,z,y);
      
  txz->SetLineColor(color);
  tyz->SetLineColor(color);
  txy->SetLineColor(color);
  
  fGraphList.push_back(txz);
  fGraphList.push_back(tyz);	
  fGraphList.push_back(txy);   
  
  XZ->cd();
  txz->Draw("same C");
  YZ->cd();
  tyz->Draw("same C");
  XY->cd();
  txy->Draw("same C");  

}


//**************************************************************
void IRunReconGlobalEvDisplay::DrawVertices(int npoints, double *x, double *y, double *z, int vertexColor, int marker) {
//**************************************************************

  if (npoints ==0) return;
  
  TGraph *txy = new TGraph(npoints,x,y);
  TGraph *txz = new TGraph(npoints,z,x);
  TGraph *tyz = new TGraph(npoints,z,y);
  
  
  txz->SetMarkerStyle(marker);
  tyz->SetMarkerStyle(marker);
  txy->SetMarkerStyle(marker);
  
  txz->SetMarkerSize(2);
  tyz->SetMarkerSize(2);
  txy->SetMarkerSize(2);

  txz->SetLineWidth(2);
  tyz->SetLineWidth(2);
  txy->SetLineWidth(2);
  
  txz->SetMarkerColor(vertexColor);
  tyz->SetMarkerColor(vertexColor);
  txy->SetMarkerColor(vertexColor);
  
  fGraphList.push_back(txz);
  fGraphList.push_back(tyz);	
  fGraphList.push_back(txy);   
  
  XZ->cd();
  txz->Draw("P");
  YZ->cd();
  tyz->Draw("P");
  XY->cd();
  txy->Draw("P");  

}


//**************************************************************
COMET::IHandle<COMET::IHitSelection>  IRunReconGlobalEvDisplay::GetTPCHits(COMET::ICOMETEvent& event) {
//**************************************************************

  //----------- Get the tpc hits. --------------------
  COMET::IHandle<COMET::IHitSelection> tpc(event.GetHitSelection("tpc"));

  if( !tpc ) { 
    std::cout << " NO TPC hits found --> Try to calibrate ..." << std::endl; 
    COMET::ICalibManager::Get().Calibrate("tpc");
  }
  tpc = event.GetHitSelection("tpc");
  if( !tpc ) { 
    std::cout << " NO TPC hits found after calibration. Skip event !!!" << std::endl; 
    return tpc;
  }

  //  std::cout << " Number of TPC hits " << tpc->size() << std::endl;
  
  return tpc;

}

//**************************************************************
COMET::IHandle<COMET::IHitSelection>  IRunReconGlobalEvDisplay::GetFGDHits(COMET::ICOMETEvent& event) {
//**************************************************************

  //------------ Get the fgd hits ---------------------
  
  COMET::IHandle<COMET::IHitSelection> fgd(event.GetHitSelection("fgd"));
  if( !fgd ) { 
    COMET::ICalibManager::Get().Calibrate("fgd");
    std::cout << " NO FGD hits found --> Try to calibrate ..." << std::endl; 
  }
  fgd = event.GetHitSelection("fgd");
  if( !fgd ) { 
    std::cout << " NO FGD hits found after calibration" << std::endl; 
    return fgd;
  }
  
  //  std::cout << " Number of FGD hits " << fgd->size() << std::endl;
  
  return fgd;
  
}

//**************************************************************
COMET::IHandle<COMET::IHitSelection> IRunReconGlobalEvDisplay::GetHitsFromComboHits(COMET::IHandle<COMET::IHitSelection> tpc) {
//**************************************************************

  if (!tpc) return COMET::IHandle<COMET::IHitSelection>();

  COMET::IHandle<COMET::IHitSelection> hits(new COMET::IHitSelection());

  for(COMET::IHitSelection::const_iterator hsel = tpc->begin(); hsel != tpc->end(); hsel++){
    COMET::IHandle<COMET::IComboHit> combo(*hsel);
    if (combo){
      for (COMET::IHitSelection::const_iterator cHit = combo->GetHits().begin();
     cHit != combo->GetHits().end();cHit++){
  if (*cHit)
    hits->push_back(*cHit);
      }
    }
  }

  return hits;
  
}


//**************************************************************
void IRunReconGlobalEvDisplay::InitExplanation() {
//**************************************************************

  double size=0.02;

  std::string s[5]={"True trajectory","Fit","Hits in fitted object","Hits in non-fitted object","Unused hits"};
  int color[5];
  color[0] = fTruthColor;
  color[1] = fObjectNodesColor;
  color[2] = fObjectHitsColor;
  color[3] = fOtherHitsColor;
  color[4] = fUnusedHitsColor;

  double deltay = 150;
  double deltax = 1700;
  double ypos_ini = 3200;
  double xpos_ini = 1000;
  double pos= ypos_ini;
  for (int i=0;i<5;i++){
    tline[i] = new TLine(xpos_ini,pos,xpos_ini+200,pos);  
    tline[i]->SetLineColor(color[i]); 
    text[i] = new TText(xpos_ini+250,pos-50,s[i].c_str());  
    text[i]->SetTextSize(size); 
    text[i]->SetTextColor(color[i]); 
    pos -=deltay;
  }

  double x[1];
  double y[1];

  x[0]=xpos_ini+deltax; y[0]=ypos_ini;
  tmarker[0] = new TGraph(1,x,y);
  tmarker[0]->SetMarkerColor(fTrueVertexColor); 
  tmarker[0]->SetMarkerSize(1.2);
  tmarker[0]->SetMarkerStyle(kFullStar);

  x[0]=xpos_ini+deltax; y[0]=ypos_ini-deltay;
  tmarker[1] = new TGraph(1,x,y);
  tmarker[1]->SetMarkerColor(fRecVertexColor); 
  tmarker[1]->SetMarkerSize(1.2);
  tmarker[1]->SetLineWidth(2);
  tmarker[1]->SetMarkerStyle(kOpenCircle);

  text[5] = new TText(xpos_ini+deltax+100,ypos_ini-50,"True primary vertex");  
  text[5]->SetTextSize(size); 
  text[5]->SetTextColor(fTrueVertexColor); 

  text[6] = new TText(xpos_ini+deltax+100,ypos_ini-deltay-50,"Rec vertex");  
  text[6]->SetTextSize(size); 
  text[6]->SetTextColor(fRecVertexColor); 

}

//**************************************************************
void IRunReconGlobalEvDisplay::DrawExplanation() {
//**************************************************************

  YZ->cd();

  for (int i=0;i<5;i++){
    tline[i]->Draw("same"); 
    text[i]->Draw("same"); 
  }


  tmarker[0]->Draw("P");
  tmarker[1]->Draw("P");

  text[5]->Draw("same"); 
  text[6]->Draw("same"); 

}

//**************************************************************
bool IRunReconGlobalEvDisplay::FindEvent(int requestedEvent, COMET::ICOMETEvent* event0, COMET::ICOMETEvent*& event) {
//**************************************************************


  int evtID = event0->GetEventId();
  int evtID2=0;

  std::cout << "evtID: " << evtID << "  " << requestedEvent << std::endl;

  COMET::ICOMETEvent* event2 = event0; 

  int readEvent = requestedEvent;
  if (evtID>=0 && evtID!=requestedEvent){
    while (evtID2 != requestedEvent){
      event2 = fInputFile->ReadEvent(readEvent);
      evtID2 = event2->GetEventId();
      std::cout << readEvent << " " << evtID2 << std::endl;
      if (evtID2 <requestedEvent)
  readEvent = readEvent + (requestedEvent-evtID2);
      else 
  readEvent = readEvent + (requestedEvent-evtID2);
    }
    //    return true;
  }


  event = event2;

  return true;

}

//**************************************************************
bool IRunReconGlobalEvDisplay::operator () (COMET::ICOMETEvent& event0) {
//**************************************************************

  // Don't do anything till we have a valid event.
  if((unsigned int)event0.GetContext().GetEvent() == 0xdeadbeef || 
     (unsigned int)event0.GetContext().GetTimeStamp() == 0xdeadbeef){    
    return true;
  }
  

  if(!geomInitialized){
    Initialized();
    GeomInit();
    geomInitialized = true;
  }

  double T0;
  double T0avg = 0;
  int NT0 = 0;

  int requestedEvent = mf->GetRequestedEvent();

  COMET::ICOMETEvent* event2;
  FindEvent(requestedEvent, &event0, event2);

  COMET::ICOMETEvent& event = *event2;

  //  COMET::ICOMETEvent& event = event0;


  //  Int_t requestedEvent=fNumber->GetIntNumber();
  int evtID = event.GetEventId();

  //  COMET::IEventFolder::GetEventFolder()->SetCurrentEvent(fInputFile->FirstEvent());


  int rCode = 0;  
  char evnum[100]; sprintf(evnum,"%4.1d", evtID);
  char rcode[100]; sprintf(rcode,"%2.1d", rCode);
  TString evstr = TString("Event No.:") + TString(evnum);
  if (isMC(event))
    evstr += TString("   Reaction code:") + TString(rcode);
  
  std::cout << " ------------------------ " << std::endl;    
  std::cout << evstr.Data() << std::endl;    
  std::cout << " ------------------------ " << std::endl;    


  //  if (evtID != 5383) return true;    

  if (isMC(event))
    fTimeOffset = 0;

  //------------- Draw the detector -------------------------
  DrawDetector();

  //------------- Draw the event information -------------------------
  DrawEventInfo(event);
  
  //----------- Draw the explanation of colors ------------------  
  DrawExplanation();

  //-------------- Run reconstruction ------------------- 
  COMET::IHandle<COMET::IAlgorithmResult>  GlobalResult = event.GetFit("ReconGlobal"); 

  // check whether ReconGlobal has been already run
  if (!GlobalResult || fForceRunGlobalRecon || fForceRunSubdetectorRecon)
    GlobalResult = fGlobalRecon->Process(event); 

  //-------------- Get Containers -------------------------
  COMET::IHandle<COMET::IAlgorithmResult>  TrackerResult = event.GetFit("ReconTracker"); 
  COMET::IHandle<COMET::IAlgorithmResult>  P0DResult     = event.GetFit("ReconP0D"); 
  COMET::IHandle<COMET::IAlgorithmResult>  ECALResult    = event.GetFit("ReconECal"); 
  COMET::IHandle<COMET::IAlgorithmResult>  P0DECALResult = event.GetFit("P0DECALPid"); 
  COMET::IHandle<COMET::IAlgorithmResult>  SMRDResult    = event.GetFit("ReconSMRD"); 


  COMET::IHandle<COMET::IHitSelection> TPCused;  
  COMET::IHandle<COMET::IHitSelection> TPCunused;
  COMET::IHandle<COMET::IHitSelection> FGDused;  
  COMET::IHandle<COMET::IHitSelection> FGDunused;
  COMET::IHandle<COMET::IHitSelection> P0Dused;  
  COMET::IHandle<COMET::IHitSelection> P0Dunused;
  COMET::IHandle<COMET::IHitSelection> ECALused; 
  COMET::IHandle<COMET::IHitSelection> ECALunused;
  COMET::IHandle<COMET::IHitSelection> P0DECALused; 
  COMET::IHandle<COMET::IHitSelection> P0DECALunused;
  COMET::IHandle<COMET::IHitSelection> SMRDused;  
  COMET::IHandle<COMET::IHitSelection> SMRDunused;


  if (TrackerResult){
    TPCused          = TrackerResult->GetHitSelection("TPC:used");
    TPCunused        = TrackerResult->GetHitSelection("TPC:unused");
    FGDused          = TrackerResult->GetHitSelection("FGD:used");
    FGDunused        = TrackerResult->GetHitSelection("FGD:unused");
  }

  if (P0DResult){
    P0Dused          = P0DResult->GetHitSelection("used");
    P0Dunused        = P0DResult->GetHitSelection("unused");
  }
  if (ECALResult){
    ECALused         = ECALResult->GetHitSelection("used");
    ECALunused       = ECALResult->GetHitSelection("unused");
  }
  if (P0DECALResult){
    P0DECALused         = P0DECALResult->GetHitSelection("used");
    P0DECALunused       = P0DECALResult->GetHitSelection("unused");
  }
  if (SMRDResult){
    SMRDused         = SMRDResult->GetHitSelection("used");
    SMRDunused       = SMRDResult->GetHitSelection("unused");
  }

  std::cout<<"Summary of hits:"<<std::endl;
    if (TPCused)
      std::cout<<" - TPC hits used:      " << TPCused->size() << std::endl;
    if (TPCunused)
      std::cout<<" - TPC hits unused:    " << TPCunused->size() << std::endl;
    if (FGDused)
      std::cout<<" - FGD hits used:      " << FGDused->size() << std::endl;
    if (FGDunused)
      std::cout<<" - FGD hits unused:    " << FGDunused->size() << std::endl;
    if (ECALused)
      std::cout<<" - ECAL hits used:     " << ECALused->size() << std::endl;
    if (ECALunused)
      std::cout<<" - ECAL hits unused:   " << ECALunused->size() << std::endl;
    if (P0DECALused)
      std::cout<<" - P0DECAL hits used:     " << P0DECALused->size() << std::endl;
    if (P0DECALunused)
      std::cout<<" - P0DECAL hits unused:   " << P0DECALunused->size() << std::endl;
    if (P0Dused)
      std::cout<<" - P0D hits used:      " << P0Dused->size() << std::endl;
    if (P0Dunused)
      std::cout<<" - P0D hits unused:    " << P0Dunused->size() << std::endl;
    if (SMRDused)
      std::cout<<" - SMRD hits used:     " << SMRDused->size() << std::endl;
    if (SMRDunused)
      std::cout<<" - SMRD hits unused:   " << SMRDunused->size() << std::endl;


  //------------- Draw all hits in FGD  -------------
  COMET::IHandle<COMET::IHitSelection> fgd = GetFGDHits(event);
  //  Draw2DHits(fgd, fUnusedHitsColor);

  //----------- draw unused hits ------------------
  if (FGDunused)
    Draw2DHits(FGDunused, fUnusedHitsColor);	  	

  if (P0Dunused)
    Draw2DHits(P0Dunused, fUnusedHitsColor);	  	

  if (ECALunused)
    Draw2DHits(ECALunused, fUnusedHitsColor);	  	

  if (P0DECALunused)
    Draw2DHits(P0DECALunused, fUnusedHitsColor);	  	

  if (SMRDunused)
    Draw2DHits(SMRDunused, fUnusedHitsColor);	  	

  //----------- draw MC-true primary vertices ------------------
  if (isMC(event))
    DrawMCVertices(event,fTrueVertexColor);

  //------------- Draw fitted objects  -------------  
  if( GlobalResult ) {    

    COMET::IHandle<COMET::IReconObjectContainer> final = GlobalResult->GetResultsContainer("final");         
    COMET::IHandle<COMET::IReconObjectContainer> vertices = GlobalResult->GetResultsContainer("vertices");         
    COMET::IHandle<COMET::IReconObjectContainer> other = GlobalResult->GetResultsContainer("other");
    
    std::cout<<"Global Recon results:"<<std::endl;
    if (vertices)
      std::cout<<" - vertices:       " << vertices->size() << std::endl;
    if (final)
      std::cout<<" - final:          " << final->size() << std::endl;
    if (other)
      std::cout<<" - other:          " << other->size() << std::endl;
    

    std::cout << std::endl;
    std::cout << "********** SUMMARY OF GLOBAL RECON *****************" << std::endl;    
    ReconPrintout::DumpContainerWithHistory(*final);   
    if (other)
      ReconPrintout::DumpContainerWithHistory(*other);   
    std::cout << "****************************************************" << std::endl;
    std::cout << std::endl;

    //----------- compute average T0 ------------------
    if( final ) {      
      for( COMET::IReconObjectContainer::iterator  tr = final->begin(); tr != final->end();tr++) { 
  COMET::IHandle<COMET::IReconBase> object = *tr;	
  
  T0 = TrackingUtils::GetPosition(object).T()+ fTimeOffset;
  if (T0 !=0){
    T0avg +=T0;
    NT0++;
  }
      }
    }

    if (NT0!=0)
      T0avg /= (1.*NT0);

    std::cout << "# T0   = " << NT0 << std::endl;    
    std::cout << "Avg T0 = " << T0avg << std::endl;
    
    //----------- draw unused TPC hits ------------------
    if (TPCunused)
      DrawTPCHits(TPCunused,T0avg,fUnusedHitsColor);

    if( final ) {      
      for( COMET::IReconObjectContainer::iterator  tr = final->begin(); tr != final->end();tr++) { 
  COMET::IHandle<COMET::IReconBase> object = *tr;	
  
  T0 = TrackingUtils::GetPosition(object).T()+ fTimeOffset;
  
  // draw MC trajectories associated to fitted objects
  if (isMC(event))
    DrawMCTrajectories(object, fTruthColor);
  
  // draw TPC hits included in fitted objects
  if (object->UsesDetector(COMET::IReconBase::kTPC))
    DrawTPCHits(object->GetHits(), T0, fObjectHitsColor);	  

  // draw 2D hits included in fitted objects
  Draw2DHits(object->GetHits(), fObjectHitsColor);	  

  // draw the results of the fitted objects
  DrawNodesInObject(object,fObjectNodesColor);		  
  
      }

    }  

    // draw reconstructed vertices
    DrawRecVertices(vertices,fRecVertexColor);

    
    //------------- Draw non-fitted objects  -------------  
    if( other ) {      
      for( COMET::IReconObjectContainer::iterator  tr = other->begin(); tr < other->end();tr++) { 
  COMET::IHandle<COMET::IReconBase> object = *tr;	
  
  // draw MC trajectories associated to non-fitted objects
  if (isMC(event))
    DrawMCTrajectories(object, fTruthColor);
  

  if (object->UsesDetector(COMET::IReconBase::kTPC)){
    T0 = TrackingUtils::GetPosition(object).T() + fTimeOffset;
    DrawTPCHits(object->GetHits(), T0, fOtherHitsColor);	  
  }
  else
    Draw2DHits(object->GetHits(), fOtherHitsColor);	  	
      }
    }  
        
  }
  else {
    std::cout << " No Fit " << std::endl; 
    // draw all TPC hits
    COMET::IHandle<COMET::IHitSelection> tpc = GetTPCHits(event);
    DrawTPCHits(tpc,0,fUnusedHitsColor);
  }

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

  YZ->Update();
  XZ->Update();
  XY->Update();


  std::cout << " ------------------------ " << std::endl;    
  std::cout << "END of " << evstr.Data() << std::endl;    
  std::cout << " ------------------------ " << std::endl;    


  theApp->Run(kTRUE);
    
  // ------------------------------------------
  for(unsigned int i = 0; i < fBoxPadList.size();i++){
    delete fBoxPadList[i];
  }  
  fBoxPadList.clear();
  
  for(unsigned int i = 0; i < fGraphList.size();i++){
    delete fGraphList[i];
  }  
  fGraphList.clear();
  
  DrawDetector();
  DrawExplanation();  

  return true;
}
  

//**************************************************************
int main(int argc, char **argv) {
//**************************************************************

  IRunReconGlobalEvDisplay userCode;
  cometEventLoop(argc,argv,userCode);
 
}