#include <iostream>

#include <TROOT.h>
#include <TProof.h>
#include <TChain.h>
#include <TString.h>
#include <TFile.h>
#include <TMacro.h>
#include <TH2D.h>

#include "rcpicker.h"

#define DOING_PROOF 1

TString proofserver("proof://stat31.phys.uconn.edu:1095");
TString xrootdbase("root://nod26.phys.uconn.edu");
TString xrootdpath("/Gluex/beamline/diamonds");
TString resultspath("INVALID");
TString runperiod("cls-12-2021");

void run_cropping(const char *prefix, Int_t seq);
void run_normalize(const char *prefix, Int_t seq);

void set_proofserver(const char *str) {
   proofserver = str;
}

void set_xrootdbase(const char *str) {
   xrootdbase = str;
}

void set_rootpath(const char *str) {
   xrootdpath = str;
}

void set_resultspath(const char *str) {
   resultspath = str;
}

void run_rcfitter(const char *prefix, Int_t seq,
                  Long64_t count=9999999999999, Long64_t start=0)
{
   if (resultspath == "INVALID") {
      std::cout << "Please set resultspath for this run, then try again."
                << std::endl
                << " example: set_resultspath(\"chess-5-2013/results\");"
                << std::endl;
      return;
   }
   TString xrootdurl;
   xrootdurl.Form("%s/%s/%s/%s_%3.3d_rocking_curves.root",
                  xrootdbase.Data(),xrootdpath.Data(),
                  resultspath.Data(),prefix,seq);
   TChain *ch=new TChain("rctree");
   ch->Add(xrootdurl);
   std::cout << xrootdurl
             << " (" << ch->GetEntries() << " steps)"
             << std::endl;
#if DOING_PROOF
   if (gProof == 0) {
      TProof::Open(proofserver);
      if (gProof == 0) {
         std::cout << "Unable to open communication with the PROOF server "
                   << std::endl
                   << proofserver
                   << std::endl;
         return;
      }
      gProof->SetParameter("PROOF_LookupOpt", "none");
      //gProof->SetParameter("PROOF_Packetizer", "TPacketizer");
      gProof->SetParameter("PROOF_MaxSlavesPerNode", (Int_t)99);
   }
   gProof->Load("Map2D.cc+O");
   gProof->Load("rcfitter.C+O");
   gProof->Load("rcpicker.C+O");
   gProof->Load("xrdport.C+O");
   gProof->Exec("xrdport_update();");
   gProof->Exec("xrdport_update();", "0");
   ch->SetProof();
#else
   gROOT->ProcessLine(".L Map2D.cc+O");
   gROOT->ProcessLine(".L rcfitter.C+O");
   gROOT->ProcessLine(".L rcpicker.C+O");
#endif
   ch->Process("rcfitter","",count,start);
#if DOING_PROOF
   TFile *foutput = new TFile;
   foutput->Open("rcfitter_output.root");
#endif
   TH2D *hbase = (TH2D*)gROOT->FindObject("hbase");
   TH2D *hamp = (TH2D*)gROOT->FindObject("hamp");
   TH2D *hmu = (TH2D*)gROOT->FindObject("hmu");
   TH2D *hsigma = (TH2D*)gROOT->FindObject("hsigma");
   TH2D *hmean = (TH2D*)gROOT->FindObject("hmean");
   TH2D *hrms = (TH2D*)gROOT->FindObject("hrms");
   TH2D *hmax = (TH2D*)gROOT->FindObject("hmax");
   TH2D *hpeak = (TH2D*)gROOT->FindObject("hpeak");
   TString resultsFile;
   resultsFile.Form("%s_%3.3d_results.root",prefix,seq);
   TFile *rFile = new TFile(resultsFile,"RECREATE");
   if (hbase) hbase->Write();
   if (hamp) hamp->Write();
   if (hmu) hmu->Write();
   if (hsigma) hsigma->Write();
   if (hmean) hmean->Write();
   if (hrms) hrms->Write();
   if (hmax) hmax->Write();
   if (hpeak) hpeak->Write();
   rFile->Close();
   delete rFile;
#if DOING_PROOF
   gProof->Close("s");
   delete gProof;
   gProof = 0;
   foutput->Close();
   delete foutput;
#endif
   delete ch;
   run_cropping(prefix, seq);
   run_normalize(prefix, seq);
}

void run_cropping(const char *prefix, Int_t seq)
{
   TString xrootdurl;
   xrootdurl.Form("%s/%s/%s/%s_%3.3d_rocking_curves.root",
                  xrootdbase.Data(),xrootdpath.Data(),
                  resultspath.Data(),prefix,seq);
   TString logfile;
   logfile.Form("%s/%s/%s",runperiod.Data(),prefix,prefix);
   rcpicker rcp(prefix, seq, xrootdurl, logfile);
   TString resultsFile;
   resultsFile.Form("%s_%3.3d_results.root",prefix,seq);
   std::cout << resultsFile << std::endl;
   TFile *rFile = new TFile(resultsFile,"UPDATE");
   TH2D *hbase  = (TH2D*)rFile->Get("hbase");
   TH2D *hamp   = (TH2D*)rFile->Get("hamp");
   TH2D *hmu    = (TH2D*)rFile->Get("hmu");
   TH2D *hsigma = (TH2D*)rFile->Get("hsigma");
   TH2D *hmean  = (TH2D*)rFile->Get("hmean");
   TH2D *hrms   = (TH2D*)rFile->Get("hrms");
   TH2D *hmax   = (TH2D*)rFile->Get("hmax");
   TH2D *hpeak  = (TH2D*)rFile->Get("hpeak");
   int ixmin=0, ixmax=0;
   int iymin=0, iymax=0;
   if (hmax) {
      rcp.Crop(hmax, &ixmin, &ixmax, &iymin, &iymax)->Write();
      std::cout << "updated hmax with ixmin,ixmax,iymin,iymax=" 
                << ixmin << "," << ixmax << "," << iymin << "," << iymax
                << std::endl;
   }
   if (hmean)  rcp.Crop(hmean, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hrms)   rcp.Crop(hrms, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hpeak)  rcp.Crop(hpeak, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hamp)   rcp.Crop(hamp, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hmu)    rcp.Crop(hmu, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hsigma) rcp.Crop(hsigma, &ixmin, &ixmax, &iymin, &iymax)->Write();
   if (hbase)  rcp.Crop(hbase, &ixmin, &ixmax, &iymin, &iymax)->Write();
   rFile->Close();
}

void run_normalize(const char *prefix, Int_t seq)
{
   TString xrootdurl;
   xrootdurl.Form("%s/%s/%s/%s_%3.3d_rocking_curves.root",
                  xrootdbase.Data(),xrootdpath.Data(),
                  resultspath.Data(),prefix,seq);
   TString logfile;
   logfile.Form("%s/%s/%s",runperiod.Data(),prefix,prefix);
   rcpicker rcp(prefix, seq, xrootdurl, logfile);
   TString resultsFile;
   resultsFile.Form("%s_%3.3d_results.root",prefix,seq);
   TFile *rFile = new TFile(resultsFile,"UPDATE");
   TH2D *hbase  = (TH2D*)rFile->Get("hbase");
   TH2D *hamp   = (TH2D*)rFile->Get("hamp");
   TH2D *hmu    = (TH2D*)rFile->Get("hmu");
   TH2D *hsigma = (TH2D*)rFile->Get("hsigma");
   TH2D *hmean  = (TH2D*)rFile->Get("hmean");
   TH2D *hrms   = (TH2D*)rFile->Get("hrms");
   TH2D *hmax   = (TH2D*)rFile->Get("hmax");
   TH2D *hpeak  = (TH2D*)rFile->Get("hpeak");
   rcp.Normalize(hbase, hamp, hmu, hsigma, hmean, hrms, hmax, hpeak);
   TH1D *rcurve  = (TH1D*)gROOT->FindObject("rcurve");
   hbase->Write();
   hamp->Write();
   hmu->Write();
   hsigma->Write();
   hmean->Write();
   hrms->Write();
   hmax->Write();
   hpeak->Write();
   rcurve->Write();
   rFile->Close();
}