//_____________________________________________________________________________
// Plug strip d-bank analysis module: below are the examples of what 
// one can do in the interactive session
//_____________________________________________________________________________

/*

  ****** to run: 
  (note naming convention: class_name = "T"+module_name+"Module"), 
  the source code of the class is stored in 

  .X cluster_ana.C
  .L TClusterAnaModule.cc+
  TStnAna x("stn.root");
  x.AddModule(new TClusterAnaModule);
  x.Run();

  ****** to display a histogram

  TClusterAnaModule* m  = (TClusterAnaModule*) x.GetModule("ClusterAna");
  m->SavePlots();


  //TCesAnaModule::Hist_t* c;
  //c = m->GetHist();
  //c->fAdcCounts->Draw();

  ****** to process one event (#11 in the ntuple) and to print the data

  x.ProcessEvent(11)
  TClusterAnaModule* m = (TClusterAnaModule*) x.GetModule("ClusterAna");
  TStnClusterBlock*   clu = m->GetClusterBlock();
  clu->Print();

*/

#include <iostream>
#include "TF1.h"
#include "TCanvas.h"
#include "Stntuple/obj/TStnHeaderBlock.hh"
#include "Stntuple/obj/TStnTrackBlock.hh"

#include "TClusterAnaModule.hh"

// ClassImp(TClusterAnaModule)
//_____________________________________________________________________________
TClusterAnaModule::TClusterAnaModule(const char* name, const char* title):
  TStnModule(name,title)
{
}

//_____________________________________________________________________________
TClusterAnaModule::~TClusterAnaModule() {
}

//_____________________________________________________________________________
void TClusterAnaModule::BookHistograms() {
  char name [200];
  char title[200];

  Delete("hist");
				// book histograms
  // ces adcs

  sprintf(name, "%s_ces_wire_adc",GetName());
  sprintf(title,"%s: Cluster: Ces Wire ADC",GetName());
  fHist.fCesWireAdc= new TH1F(name, title,500 ,0.001,500);
  AddHistogram(fHist.fCesWireAdc);

  sprintf(name, "%s_ces_strp_adc",GetName());
  sprintf(title,"%s: Cluster: Ces Strp ADC",GetName());
  fHist.fCesStrpAdc= new TH1F(name, title,500 ,0.001,500);
  AddHistogram(fHist.fCesStrpAdc);

  sprintf(name, "%s_ces_wire_avg_adc",GetName());
  sprintf(title,"%s: Cluster: Ces Wire Avg ADC",GetName());
  fHist.fCesWireAAdc= new TProfile(name, title, 128,-0.5,127.5);
  AddHistogram(fHist.fCesWireAAdc);

  sprintf(name, "%s_ces_strp_avg_adc",GetName());
  sprintf(title,"%s: Cluster: Ces Strp Avg ADC",GetName());
  fHist.fCesStrpAAdc= new TProfile(name, title,500 ,-0.5,254.5);
  AddHistogram(fHist.fCesStrpAAdc);

  sprintf(name, "%s_ces_wire_occ",GetName());
  sprintf(title,"%s: Cluster: Ces Wire Occ",GetName());
  fHist.fCesWireOcc = new TProfile(name, title,128 ,-0.5,127.5);
  AddHistogram(fHist.fCesWireOcc);

  sprintf(name, "%s_ces_strp_occ",GetName());
  sprintf(title,"%s: Cluster: Ces Strp Occ",GetName());
  fHist.fCesStrpOcc = new TProfile(name, title,128 ,-0.5,254.5);
  AddHistogram(fHist.fCesStrpOcc);

  sprintf(name, "%s_ces_wire_gain",GetName());
  sprintf(title,"%s: Cluster: Ces Wire Gain",GetName());
  fHist.fCesWireGain = new TProfile(name, title,128 ,-0.5,127.5);
  AddHistogram(fHist.fCesWireGain);

  sprintf(name, "%s_ces_strp_gain",GetName());
  sprintf(title,"%s: Cluster: Ces Strp Gain",GetName());
  fHist.fCesStrpGain = new TProfile(name, title,128 ,-0.5,254.5);
  AddHistogram(fHist.fCesStrpGain);

  //cpr adc's

  sprintf(name, "%s_cpr_wire_adc",GetName());
  sprintf(title,"%s: Cluster: Cpr Wire ADC",GetName());
  fHist.fCprWireAdc= new TH1F(name, title,500 ,0.001,500);
  AddHistogram(fHist.fCprWireAdc);

  sprintf(name, "%s_cpr_wire_avg_adc",GetName());
  sprintf(title,"%s: Cluster: Cpr Wire Avg ADC",GetName());
  fHist.fCprWireAAdc= new TProfile(name, title,500 ,-0.5,64.5);
  AddHistogram(fHist.fCprWireAAdc);

  sprintf(name, "%s_cpr_wire_occ",GetName());
  sprintf(title,"%s: Cluster: Cpr Wire Occ",GetName());
  fHist.fCprWireOcc = new TProfile(name, title,128 ,-0.5,64.5);
  AddHistogram(fHist.fCprWireOcc);

  sprintf(name, "%s_cpr_wire_gain",GetName());
  sprintf(title,"%s: Cluster: Cpr Wire Gain",GetName());
  fHist.fCprWireGain = new TProfile(name, title,128 ,-0.5,64.5);
  AddHistogram(fHist.fCprWireGain);


  // ces clusters
  fHist.fNCesUnbWire = new 
    TH1F("Clus_nces_wire_unb","Cluster NCes Wire Unb",19,-0.5,20.5);
  AddHistogram(fHist.fNCesUnbWire);
  fHist.fCesUnbWireWed = new 
    TH1F("Clus_ces_unb_wire_wed","Cluster Ces Unb Wire Wed",48,-0.5,47.5);
  AddHistogram(fHist.fCesUnbWireWed);
  fHist.fCesUnbWireE = new 
    TH1F("Clus_ces_unb_wire_e","Cluster Ces Unb Wire E",50,0.,30.0);
  AddHistogram(fHist.fCesUnbWireE);
  fHist.fCesUnbX = new 
    TH1F("Clus_ces_unb_x","Cluster Ces Unb X",50,-30.0,30.0);
  AddHistogram(fHist.fCesUnbX);
  fHist.fCesUnbWireSigma = new 
    TH1F("Clus_ces_unb_wire_sigma","Cluster Ces Unb Wire Sigma",50,0.0,10.0);
  AddHistogram(fHist.fCesUnbWireSigma);
  fHist.fCesUnbWireChiSq = new 
    TH1F("Clus_ces_unb_wire_chisq","Cluster Ces Unb Wire ChiSq",50,0.0,30.0);
  AddHistogram(fHist.fCesUnbWireChiSq);

  fHist.fNCesUnbStrp = new 
    TH1F("Clus_nces_strp_unb","Cluster NCes Strp Unb",19,-0.5,20.5);
  AddHistogram(fHist.fNCesUnbStrp);
  fHist.fCesUnbStrpWed = new 
    TH1F("Clus_ces_unb_strp_wed","Cluster Ces Unb Strp Wed",48,-0.5,47.5);
  AddHistogram(fHist.fCesUnbStrpWed);
  fHist.fCesUnbStrpE = new 
    TH1F("Clus_ces_unb_strp_e","Cluster Ces Unb Strp E",50,0.,30.);
  AddHistogram(fHist.fCesUnbStrpE);
  fHist.fCesUnbZ = new 
    TH1F("Clus_ces_unb_z","Cluster Ces Unb Z",50,0.0,250.0);
  AddHistogram(fHist.fCesUnbZ);
  fHist.fCesUnbStrpSigma = new 
    TH1F("Clus_ces_unb_strp_sigma","Cluster Ces Unb Strp Sigma",50,0.0,10.0);
  AddHistogram(fHist.fCesUnbStrpSigma);
  fHist.fCesUnbStrpChiSq = new 
    TH1F("Clus_ces_unb_strp_chisq","Cluster Ces Unb Strp ChiSq",50,0.0,30.0);
  AddHistogram(fHist.fCesUnbStrpChiSq);


  // cpr clusters
  fHist.fNCprUnb= new TH1F("Clus_ncpr_unb","Cluster NCpr Unb",19,-0.5,20.5);
  AddHistogram(fHist.fNCprUnb);
  fHist.fCprUnbWed= new 
    TH1F("Clus_cpr_unb_wed","Cluster Cpr Unb Wedge",48,-0.5,47.5);
  AddHistogram(fHist.fCprUnbWed);
  fHist.fCprUnbE= new 
    TH1F("Clus_cpr_unb_e","Cluster Cpr Unb Q",100,0.0,45000.0);
  AddHistogram(fHist.fCprUnbE);
  fHist.fCprUnbX= new TH1F("Clus_cpr_unb_x","Cluster Cpr Unb X",50,-20.0,20.0);
  AddHistogram(fHist.fCprUnbX);
  fHist.fCprUnbWireSigma = new 
           TH1F("Clus_cpr_unb_sig","Cluster Cpr Unb Sig",19,-0.5,2000.5);
  AddHistogram(fHist.fCprUnbWireSigma);
  fHist.fCprUnbWireChiSq = new 
           TH1F("Clus_cpr_unb_chi2","Cluster Cpr Unb ChiSq",19,-0.5,20.5);
  AddHistogram(fHist.fCprUnbWireChiSq);

}


//_____________________________________________________________________________
int TClusterAnaModule::BeginJob() {

  // register the data blocks and book the histograms

  RegisterDataBlock("CesDataBlock","TCesDataBlock"   ,&fCesDataBlock);
  RegisterDataBlock("CprDataBlock","TCprDataBlock"   ,&fCprDataBlock);
  RegisterDataBlock("ClusterBlock","TStnClusterBlock",&fClusterBlock);
  
  BookHistograms();

  return 0;
}


//_____________________________________________________________________________
int TClusterAnaModule::BeginRun() {
  return 0;
}

//_____________________________________________________________________________
int TClusterAnaModule::Event(int ientry) {

  // each of these is a TBranch, move the brach to this event
  fClusterBlock->GetEntry(ientry);

  if (fCesDataBlock) {
    fCesDataBlock->GetEntry(ientry);
    for (int is=0; is<2; is++) {
      for (int iw=0; iw<24; iw++) {
	if (fCesDataBlock->GetNHitStrips(is,iw) > 0) {
	  for (int ich=0; ich<128; ich++) {
	    int adc_counts = fCesDataBlock->GetStripData(is,iw,ich);
	    fHist.fCesStrpAdc->Fill(adc_counts);
	    float x = (iw-14)*128+ich;
	    fHist.fCesStrpAAdc->Fill(x,adc_counts);
	    float w = ( adc_counts>120 ? 1 : 0);
	    fHist.fCesStrpOcc->Fill(x,w);
	    if(adc_counts>120.0) fHist.fCesStrpGain->Fill(x,adc_counts);
	  }
	}

	if (fCesDataBlock->GetNHitWires(is,iw) > 0) {
	  for (int ich=0; ich<64; ich++) {
	    int adc_counts = fCesDataBlock->GetWireData(is,iw,ich);
	    fHist.fCesWireAdc->Fill(adc_counts);
	    float x = (iw-14)*64+ich;
	    fHist.fCesWireAAdc->Fill(x,adc_counts);
	    float w = ( adc_counts>120 ? 1 : 0);
	    fHist.fCesWireOcc->Fill(x,w);
	    if(adc_counts>120.0) fHist.fCesWireGain->Fill(x,adc_counts);
	  }
	}
      }
    }
  }
//-----------------------------------------------------------------------------
// CPR wire info
//-----------------------------------------------------------------------------
  if (fCprDataBlock) {
    fCprDataBlock->GetEntry(ientry);
    for (int is=0; is<2; is++) {
      for (int iw=0; iw<24; iw++) {
	if (fCprDataBlock->GetNHitWires(is,iw) > 0) {
	  for (int ich=0; ich<32; ich++) {
	    int adc_counts = fCprDataBlock->GetWireData(is,iw,ich);
	    fHist.fCprWireAdc->Fill(adc_counts);
	    float x = (iw-14)*32+ich;
	    fHist.fCprWireAAdc->Fill(x,adc_counts);
	    float w = ( adc_counts>120 ? 1 : 0);
	    fHist.fCprWireOcc->Fill(x,w);
	    if(adc_counts>120.0) fHist.fCprWireGain->Fill(x,adc_counts);
	  }
	}
      }
    }
  }
//-----------------------------------------------------------------------------
// clusters
//-----------------------------------------------------------------------------
  int nPerWed[48];
  for(int i=0;i<48;i++) nPerWed[i]=0;
  int n = fClusterBlock->NCesUnbClusters();
  for(int i=0;i<fClusterBlock->NCesUnbClusters(); i++) {
    TCesCluster& ces = *(fClusterBlock->CesUnbCluster(i));
    if(ces.View()==0) {
      int ind = (ces.Side()*24) + ces.Wedge();
      nPerWed[ind]++;
      fHist.fCesUnbWireWed->Fill(ind);
      fHist.fCesUnbWireE->Fill(ces.Energy());
      fHist.fCesUnbX->Fill(ces.LocalCoord());
      fHist.fCesUnbWireSigma->Fill(ces.Sigma());
      fHist.fCesUnbWireChiSq->Fill(ces.ChiSq());
    }
  }
  for(int i=0;i<48; i++) {
    fHist.fNCesUnbWire->Fill(nPerWed[i]);
  }

  for(int i=0;i<48;i++) nPerWed[i]=0;
  n = fClusterBlock->NCesUnbClusters();
  for(int i=0;i<fClusterBlock->NCesUnbClusters(); i++) {
    TCesCluster& ces = *(fClusterBlock->CesUnbCluster(i));
    if(ces.View()==1) {
      int ind = (ces.Side()*24) + ces.Wedge();
      nPerWed[ind]++;
      fHist.fCesUnbStrpWed->Fill(ind);
      fHist.fCesUnbStrpE->Fill(ces.Energy());
      fHist.fCesUnbZ->Fill(ces.LocalCoord());
      fHist.fCesUnbStrpSigma->Fill(ces.Sigma());
      fHist.fCesUnbStrpChiSq->Fill(ces.ChiSq());
    }
  }
  for(int i=0;i<48; i++) {
    fHist.fNCesUnbStrp->Fill(nPerWed[i]);
  }


  // cpr clusters
  for(int i=0;i<48;i++) nPerWed[i]=0;
  n = fClusterBlock->NCprUnbClusters();
  for(int i=0;i<n; i++) {
    TCprCluster& cpr = *(fClusterBlock->CprUnbCluster(i));
    int ind = (cpr.Side()*24) + cpr.Wedge();
    nPerWed[ind]++;
    fHist.fCprUnbWed->Fill(ind);
    fHist.fCprUnbE->Fill(cpr.Charge());
    fHist.fCprUnbX->Fill(cpr.X());
    fHist.fCprUnbWireSigma->Fill(cpr.SigmaX());
    //fHist.fCprUnbWireChiSq->Fill(cpr.());
  }
  for(int i=0;i<48; i++) {
    fHist.fNCprUnb->Fill(nPerWed[i]);
  }

  return 0;
  }


//_____________________________________________________________________________
void TClusterAnaModule::DisplayEvent() {

  /*
  char name[]  = "TClusterAnaModule_DisplayEvent_canvas";
  char title[] = "West plug Cluster events";

  TCanvas* c = (TCanvas*) gROOT->FindObject(name);

  if (! c) {
    c = new TCanvas(name, title,0,0,800,600);
    c->Divide(2,2);
  }

  for (int k=0; k<2; k++) {
    c->cd(2*k+1);
    fHist.fStripProfile[k]->Draw();
    c->cd(2*k+2);
    fHist.fWireProfile[k]->Draw();
  }
  c->Modified();
  c->Update();
  */
}

void TClusterAnaModule::SavePlots() {

  TCanvas * c1= new TCanvas("Clus_page1","Cluster Page 1",100,100,600,600);
  c1->Divide(2,3);
  c1->cd(1); fHist.fCesWireAdc->Draw();
  c1->cd(2); fHist.fCesWireAAdc->Draw();
  c1->cd(3); fHist.fCesStrpAdc->Draw();
  fHist.fCesStrpAAdc->SetMinimum(-200);
  fHist.fCesStrpAAdc->SetMaximum( 200);
  c1->cd(4); fHist.fCesStrpAAdc->Draw();
  c1->cd(5); fHist.fCprWireAdc->Draw();
  c1->cd(6); fHist.fCprWireAAdc->Draw();

  TCanvas * c2= new TCanvas("Clus_page2","Cluster Page 2",150,150,600,600);
  c2->Divide(2,3);
  c2->cd(1); fHist.fCesWireOcc->Draw();
  c2->cd(2); fHist.fCesWireGain->Draw();
  c2->cd(3); fHist.fCesStrpOcc->Draw();
  c2->cd(4); fHist.fCesStrpGain->Draw();
  c2->cd(5); fHist.fCprWireOcc->Draw();
  c2->cd(6); fHist.fCprWireGain->Draw();

  TCanvas * c3= new TCanvas("Clus_page3","Cluster Page 3",200,200,600,600);
  c3->Divide(2,3);
  c3->cd(1); fHist.fNCprUnb->Draw();
  c3->cd(2); fHist.fCprUnbWed->Draw();
  c3->cd(3); fHist.fCprUnbE->Draw();
  c3->cd(4); fHist.fCprUnbX->Draw();
  c3->cd(5); fHist.fCprUnbWireSigma->Draw();
  c3->cd(6); fHist.fCprUnbWireChiSq->Draw();

  TCanvas * c4= new TCanvas("Clus_page4","Cluster Page 4",250,250,600,600);
  c4->Divide(2,3);
  c4->cd(1); fHist.fNCesUnbWire->Draw();
  c4->cd(2); fHist.fNCesUnbStrp->Draw();
  c4->cd(3); fHist.fCesUnbWireWed->Draw();
  c4->cd(4); fHist.fCesUnbStrpWed->Draw();
  c4->cd(5); fHist.fCesUnbWireE->Draw();
  c4->cd(6); fHist.fCesUnbStrpE->Draw();

  TCanvas * c5= new TCanvas("Clus_page5","Cluster Page 5",300,300,600,600);
  c5->Divide(2,3);
  c5->cd(1); fHist.fCesUnbX->Draw();
  c5->cd(2); fHist.fCesUnbZ->Draw();
  c5->cd(3); fHist.fCesUnbWireSigma->Draw();
  c5->cd(4); fHist.fCesUnbStrpSigma->Draw();
  c5->cd(5); fHist.fCesUnbWireChiSq->Draw();
  c5->cd(6); fHist.fCesUnbStrpChiSq->Draw();

}

//_____________________________________________________________________________
int TClusterAnaModule::EndJob() {
  printf("----- end job: ---- %s\n",GetName());
  return 0;
}