// created    : May 23rd, 2003
// by         : Gero Flucke
// last update: $Date: 2005/12/29 16:02:26 $
// by         : $Author: flucke $

#include "TObjArray.h"
#include "TObjString.h"
#include "TH1.h"
#include "TF1.h"
#include <TGraphAsymmErrors.h>
#include "TLine.h"
#include "TCanvas.h"
#include <TLegend.h>

#include <iostream>

#if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
#include "H1Steering/H1Steer.h"
#endif
#include "GFAnalysis/GFDstarCompare.h"
#include "GFAnalysis/GFDstarHistsAnalysisData.h"
#include "GFAnalysis/GFDstarHistsAnalysisMc.h"
#include "GFAnalysis/GFMultiDmFitter.h"

#include "GFUtils/GFHistManip.h"
#include "GFUtils/GFFuncs.h"
#include "GFUtils/GFHistArray.h"
#include "GFUtils/GFHistManager.h"
#include "GFUtils/GFDstarRunPeriods.h"

#include "/afs/desy.de/user/f/flucke/h1/root/DstarDmFitter/DstarDmFitter.h"

using namespace std;

typedef vector<GFDstarHistsAnalysisMc*>::iterator McIter;
typedef vector<GFDstarHistsAnalysisMc*>::const_iterator ConstMcIter;

typedef vector<GFDstarHistsAnalysisData*>::iterator DataIter;
typedef vector<GFDstarHistsAnalysisData*>::const_iterator ConstDataIter;

ClassImp(GFDstarCompare)

 GFDstarCompare::GFDstarCompare()
{
  fHistManager = new GFHistManager;
}

 GFDstarCompare::GFDstarCompare(const char* filePrefix)
{
  fHistManager = new GFHistManager;

  //  TString dir("/data/usr/flucke/dstar/data/oo/histFiles/");
//   TString name(filePrefix);
//   name += "DstarAnalysisHists";

  GFDstarHistsAnalysisData *newData = new GFDstarHistsAnalysisData(filePrefix, "All");
  newData->SetBatch();
//   newData->SetPeriod("Data (fit)");
//   newData->GetMultiDmFitter()->SetMergeNBins(2);
  fDatas.push_back(newData);

//   newData = new GFDstarHistsAnalysisData(filePrefix, "All");
//   newData->SetBatch();
//   newData->SetPeriod("Data (SR)");
//   newData->GetMultiDmFitter()->SetFitOrNotFlag(2);
//   fDatas.push_back(newData);

  GFDstarHistsAnalysisMc *newMc = NULL;

//   newMc = new GFDstarHistsAnalysisMcMore(filePrefix);// before 'per lumi'
  newMc = new GFDstarHistsAnalysisMc(filePrefix, "Pythia");
  newMc->SetBatch();
  fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharm.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharm99ePlus.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmA.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmB.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmC.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmD.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmAB.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "ResCharm.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "ResCharm99ePlus.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "ResCharmA.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "ResCharmB.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CharmEx.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CharmEx99ePlus.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CharmExA.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CharmExB.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharm.root", 
// 					 dir + name + "ResCharmA.root", 
// 					 dir + name + "CharmEx.root");
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirCharmAB.root", 
// 					 dir + name + "ResCharmA.root", 
// 					 dir + name + "CharmEx.root");
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "DirBeauty.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "ResBeauty.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

// before 'per lumi'
//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "Cascade.root", NULL, NULL);

  newMc = new GFDstarHistsAnalysisMc(filePrefix, "Cascade");
  newMc->SetBatch();
  fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore("GFData/DstarAnalysisHistsCascade.root", NULL, NULL);
//   newMc->SetPeriod("CascadeWeighted");
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "Cascade99ePlus.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CascadeA.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CascadeB.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   newMc = new GFDstarHistsAnalysisMcMore(dir + name + "CascadeC.root", NULL, NULL);
//   newMc->SetBatch();
//   fMonteCarlo.push_back(newMc);

//   for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
//     (*i)->GetMultiDmFitter()->SetFitOrNotFlag(2);
//   }
}

 GFDstarCompare::GFDstarCompare(const char* data, const char* mc)
{
  // take full filenames, different MC may be separated by ','
  fHistManager = new GFHistManager;

  if(data){
    GFDstarHistsAnalysisData *newData = new GFDstarHistsAnalysisData(data);
    newData->SetBatch();
    fDatas.push_back(newData);
  }

#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,0)
  TObjArray *mcArray = TString(mc).Tokenize(",");// mcArray is owner!
#else
  TObjArray *mcArray = new TObjArray;
  mcArray->SetOwner();
  H1Steer help;
  help.SplitString(mc, ",", mcArray);
#endif
  for(Int_t i = 0; i < mcArray->GetEntriesFast(); ++i){
    GFDstarHistsAnalysisMc *newMc = new GFDstarHistsAnalysisMc(mcArray->At(i)->GetName());
    newMc->SetBatch();
    fMonteCarlo.push_back(newMc);
  }
  delete mcArray;
}


 GFDstarCompare::~GFDstarCompare()
{
  delete fHistManager;

  for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
    delete *i;
  }

  for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
    delete *i;
  }

}

 void GFDstarCompare::AdoptData(GFDstarHistsAnalysisData *data)
{
  if(data) {
    data->SetBatch();
    fDatas.push_back(data);
  }
}

 void GFDstarCompare::AdoptMc(GFDstarHistsAnalysisMc *mc)
{
  if(mc) {
    mc->SetBatch();
    fMonteCarlo.push_back(mc);
  }
}


 TH1* GFDstarCompare::L4Eff(const char *var, Int_t dstarL4)
{
  // return first data hist
  TH1* firstDataHist = NULL;
  TObjArray vars;
  if(var) {
    vars.Add(new TObjString(var));
  } else {
    vars.Add(new TObjString("pt"));
    vars.Add(new TObjString("eta"));
    vars.Add(new TObjString("phi"));
    vars.Add(new TObjString("wGammaP"));
    vars.Add(new TObjString("y"));
    vars.Add(new TObjString("ptSumEt"));
    vars.Add(new TObjString("nTrack"));
  }
  const Int_t layer = 0;

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Long_t numData = 0;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      ++numData;
      TH1* hist = (*i)->CreateL4Eff(vars[iVar]->GetName(), dstarL4);
      TGraphAsymmErrors* graph = (*i)->TmpGraph();
      TString label = (*i)->GetPeriod();
      if(numData == 1) {
	fHistManager->AddHist(hist, layer, label, "LP");
	firstDataHist = hist;
      } else fHistManager->AddHistSame(hist, layer, iVar, label, "LP");
      if(graph){
	for(Int_t i2 = 1; i2 <= hist->GetNbinsX(); ++i2) hist->SetBinError(i2, 1.e-7);
	GFHistManager::MakeDifferentStyle(fHistManager->GetHistsOf(layer, iVar));
	hist->TAttMarker::Copy(*graph);
	hist->TAttLine::Copy(*graph);
	fHistManager->AddObject(graph, layer, iVar, "P"); 
      }
    }

    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateL4Eff(vars[iVar]->GetName(), dstarL4);
      TGraphAsymmErrors* graph = (*i)->TmpGraph();

      if(numData == 0) fHistManager->AddHist(hist, layer, (*i)->GetPeriod());
      else fHistManager->AddHistSame(hist, layer, iVar, (*i)->GetPeriod());
      if(graph){
	for(Int_t i2 = 1; i2 <= hist->GetNbinsX(); ++i2) hist->SetBinError(i2, 1.e-7);
	GFHistManager::MakeDifferentStyle(fHistManager->GetHistsOf(layer, iVar));
	hist->TAttMarker::Copy(*graph);
	hist->TAttLine::Copy(*graph);
	fHistManager->AddObject(graph, layer, iVar, "P"); 
      }
    }
  }

  fHistManager->SetHistsMinMax(0.75, 1.04);
  fHistManager->Draw();
  vars.Delete();

  return firstDataHist;
}

//_______________________________________________________________
 void GFDstarCompare::RecEff(const char* var)
{
  TObjArray vars;
  if(var){
    vars.Add(new TObjString(var));
  } else {
    vars.Add(new TObjString("pt"));
    vars.Add(new TObjString("eta"));
    vars.Add(new TObjString("phi"));
    vars.Add(new TObjString("wGammaP"));
    vars.Add(new TObjString("y"));
    vars.Add(new TObjString("ptSumEt"));
    vars.Add(new TObjString("nTrack"));
  }

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t numMc = 0;
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateRecEff(vars[iVar]->GetName());
      if(numMc == 0){
	fHistManager->AddHist(hist, 0, (*i)->GetPeriod());
      } else {
	fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
      }
      ++numMc;
    }
  }

  fHistManager->Draw();
  vars.Delete();
}

//_______________________________________________________________
 void GFDstarCompare::DetAcceptance(const char* var)
{
  TObjArray vars;
  if(var){
    vars.Add(new TObjString(var));
  } else {
    vars.Add(new TObjString("pt"));
    vars.Add(new TObjString("eta"));
    vars.Add(new TObjString("phi"));
    vars.Add(new TObjString("wGammaP"));
    vars.Add(new TObjString("y"));
    vars.Add(new TObjString("ptSumEt"));
    vars.Add(new TObjString("nTrack"));
  }

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t numMc = 0;
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateDetAcceptance(vars[iVar]->GetName());
      if(numMc == 0){
	fHistManager->AddHist(hist, 0, (*i)->GetPeriod());
      } else {
	fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
      }
      TGraph* gr = (*i)->TmpGraph();
      fHistManager->AddObject(gr, 0, iVar, "P");
      ++numMc;
    }
  }

  fHistManager->Draw();
  for(Int_t l = 0; l < fHistManager->GetNumLayers(); ++l){
    for(Int_t nH = 0; nH < fHistManager->GetNumHistsOf(l); ++nH){
      GFHistArray *hists = fHistManager->GetHistsOf(l, nH);
      TIter graphs(fHistManager->GetObjectsOf(l, nH));
      Int_t i = 0;
      while(TGraph *gr = static_cast<TGraph*>(graphs())){
	hists->At(i)->TAttLine::Copy(*gr);
	hists->At(i)->TAttMarker::Copy(*gr);
	++i;
      }
      
    }
  }
  fHistManager->Update();
  vars.Delete();
}

//_______________________________________________________________
 void GFDstarCompare::HadCorr(const char *var, const char *forwBack)
{
  TObjArray vars;
  if(var){
    vars.Add(new TObjString(var));
  } else {
    vars.Add(new TObjString("PtDs"));
    vars.Add(new TObjString("PtJet"));
    vars.Add(new TObjString("EtaDs"));
    vars.Add(new TObjString("EtaJet"));
    vars.Add(new TObjString("Deta"));
    vars.Add(new TObjString("Dphi"));
    vars.Add(new TObjString("PtDsJet"));
    vars.Add(new TObjString("MDsJet"));
    vars.Add(new TObjString("xGam"));
  }

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t numMc = 0;
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateHadCorr(vars[iVar]->GetName(), forwBack);
      if(numMc == 0){
	fHistManager->AddHist(hist, 0, (*i)->GetPeriod());
      } else {
	fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
      }
//       TGraph* gr = (*i)->TmpGraph();
//       if(gr) fHistManager->AddObject(gr, 0, iVar, "P");
      ++numMc;
    }
  }

  fHistManager->Draw();
//   for(Int_t l = 0; l < fHistManager->GetNumLayers(); ++l){
//     for(Int_t nH = 0; nH < fHistManager->GetNumHistsOf(l); ++nH){
//       GFHistArray *hists = fHistManager->GetHistsOf(l, nH);
//       TIter graphs(fHistManager->GetObjectsOf(l, nH));
//       Int_t i = 0;
//       while(TGraph *gr = static_cast<TGraph*>(graphs())){
// 	if(!gr->InheritsFrom(TGraph::Class())) continue;
// 	hists->At(i)->TAttLine::Copy(*gr);
// 	hists->At(i)->TAttMarker::Copy(*gr);
// 	++i;
//       }
//     }
//   }
  fHistManager->Update();
  vars.Delete();
}

//_______________________________________________________________
 void GFDstarCompare::Accept(const char* var)
{
  TObjArray vars;
  if(var){
    vars.Add(new TObjString(var));
  } else {
    vars.Add(new TObjString("pt"));
    vars.Add(new TObjString("eta"));
    vars.Add(new TObjString("phi"));
    vars.Add(new TObjString("wGammaP"));
    vars.Add(new TObjString("y"));
    vars.Add(new TObjString("ptSumEt"));
    vars.Add(new TObjString("nTrack"));
  }

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t numMc = 0;
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateAccept(vars[iVar]->GetName());
      if(numMc == 0){
	fHistManager->AddHist(hist, 0, (*i)->GetPeriod());
      } else {
	fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
      }
      ++numMc;
    }
  }

  fHistManager->Draw();
  vars.Delete();
}

//_______________________________________________________________
 void GFDstarCompare::DstarTracks()
{
  fHistManager->Clear();
  TObjArray vars;
  vars.Add(new TObjString("PtK"));
  vars.Add(new TObjString("ThetaK"));
  vars.Add(new TObjString("PhiK"));
  vars.Add(new TObjString("LengthK"));
  vars.Add(new TObjString("StartK"));
  vars.Add(new TObjString("PtPi"));
  vars.Add(new TObjString("ThetaPi"));
  vars.Add(new TObjString("PhiPi"));
  vars.Add(new TObjString("LengthPi"));
  vars.Add(new TObjString("StartPi"));
  vars.Add(new TObjString("PtPis"));
  vars.Add(new TObjString("ThetaPis"));
  vars.Add(new TObjString("PhiPis"));
  vars.Add(new TObjString("LengthPis"));
  vars.Add(new TObjString("StartPis"));

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Long_t numData = 0;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      TH1* hist = (*i)->DrawTrack(vars[iVar]->GetName());
      TString label = (*i)->GetPeriod(); 
      if(numData == 0) fHistManager->AddHist(hist, 0, label);
      else fHistManager->AddHistSame(hist, 0, iVar, (label += "_") +=numData);
      ++numData;
    }

    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->DrawTrack(vars[iVar]->GetName());
      fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
    }
  }
  vars.Delete();
  fHistManager->Draw();
}

//_______________________________________________________________
 void GFDstarCompare::TrigEff(Int_t TE)
{
  fHistManager->Clear();
  TObjArray vars;
  vars.Add(new TObjString("pt"));
  vars.Add(new TObjString("eta"));
  vars.Add(new TObjString("phi"));
  vars.Add(new TObjString("wGammaP"));
  //   vars.Add(new TObjString("y")); only data...
  // FIXME: more hists!

  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Long_t numData = 0;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      TH1* hist = (*i)->CreateTrigEff(vars[iVar]->GetName(), TE);
      TString label = (*i)->GetPeriod(); 
      if(numData == 0) fHistManager->AddHist(hist, 0, label);
      else fHistManager->AddHistSame(hist, 0, iVar, (label += "_") +=numData);
      ++numData;
    }

    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->CreateTrigEff(vars[iVar]->GetName(),TE);
      fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
    }
  }
  vars.Delete();
  fHistManager->Draw();
}

//_______________________________________________________________
 void GFDstarCompare::JetProfile(Int_t flagBins, Int_t etaPhi, Bool_t second)
{
  // flagBins = 0: overall only, 1: overall and in pt/eta bins, 2: only pt/eta bins
  // if(second) draw also hists for 2nd jet
  // etaPhi < 0: phi profiles, > 0 eta profiles, == 0 both

  TObjArray vars;
  Int_t nOverall = 0;
  if(flagBins < 2){
    if(etaPhi >= 0){
      vars.Add(new TObjString("jet1Eta"));
      ++nOverall;
    }
    if(etaPhi <= 0){
      vars.Add(new TObjString("jet1Phi"));
      ++nOverall;
    }
    if(second){// here ignore etaPhi...
      vars.Add(new TObjString("jet2Eta"));
      vars.Add(new TObjString("jet2Phi"));
      nOverall += 2;
    }
  }

  
  Int_t nPt = 0;
  Int_t nEta = 0;
  if(flagBins > 0){
    TH1* histJetPt = fDatas[0]->GetHist("dstar1JetPtJet", "dstarJets");
    const Int_t nPtBins = (histJetPt ? histJetPt->GetNbinsY(): 0);
    for(Int_t i = 1; i <= nPtBins; ++i){
      if(etaPhi >= 0) {
	vars.Add(new TObjString(Form("jet1EtaPt%d",i)));
	++nPt;
      }
      if(etaPhi <= 0) {
	vars.Add(new TObjString(Form("jet1PhiPt%d",i)));
	++nPt;
      }
      if(second){ // ignore etaPhi
	vars.Add(new TObjString(Form("jet2EtaPt%d",i)));
	vars.Add(new TObjString(Form("jet2PhiPt%d",i)));
	nPt += 2;
      }
    }
    TH1* histJetEta = fDatas[0]->GetHist("dstar1JetEtaJet", "dstarJets");
    const Int_t nEtaBins = (histJetEta ? histJetEta->GetNbinsY(): 0);
    for(Int_t i = 1; i <= nEtaBins; ++i){
      if(etaPhi >= 0){
	vars.Add(new TObjString(Form("jet1EtaEta%d",i)));
	++nEta;
      }
      if(etaPhi <= 0){
	vars.Add(new TObjString(Form("jet1PhiEta%d",i)));
	++nEta;
      }
      if(second){ // ignore etaPhi
	vars.Add(new TObjString(Form("jet2EtaEta%d",i)));
	vars.Add(new TObjString(Form("jet2PhiEta%d",i)));
	nEta += 2;
      }
    }
    }
  Int_t nXgamma = 0;
  if(etaPhi >= 0) {vars.Add(new TObjString("jet1EtaxGam1")); ++nXgamma;}
  if(etaPhi <= 0) {vars.Add(new TObjString("jet1PhixGam1")); ++nXgamma;}
  if(second){ // ignore etaPhi
    vars.Add(new TObjString("jet2EtaxGam1"));
    vars.Add(new TObjString("jet2PhixGam1"));
    nXgamma += 2;
  }
  if(etaPhi >= 0) {vars.Add(new TObjString("jet1EtaxGam2")); ++nXgamma;}
  if(etaPhi <= 0) {vars.Add(new TObjString("jet1PhixGam2")); ++nXgamma;}
  if(second){ // ignore etaPhi
    vars.Add(new TObjString("jet2EtaxGam2"));
    vars.Add(new TObjString("jet2PhixGam2"));
    nXgamma += 2;
  }

  const Int_t mcFitFlag = 2;
  const char* nonDataOpt = "M";
  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t nLayer = 0;
    Int_t nPad = iVar;
    if(iVar >= nOverall){
      nLayer += (nOverall ? 1 : 0);
      nPad = iVar - nOverall;
      if(iVar >= nOverall + nPt) {
	nLayer += (nPt ? 1 : 0);
	nPad = iVar - nOverall - nPt;
	if(iVar >= nOverall + nPt + nEta) {
	  nLayer += (nEta ? 1 : 0);
	  nPad = iVar - nOverall - nPt - nEta;
	}
      }
    }
    Long_t numData = 0;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      TString dataFitOpt = (*i)->GetMultiDmFitter()->GetDstarDmFitter()->GetFitOption();
      Bool_t replaced = kFALSE;
      if(dataFitOpt.Contains(nonDataOpt)) {
	dataFitOpt.ReplaceAll(nonDataOpt, NULL);
	(*i)->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(dataFitOpt);
	replaced = kTRUE;
      }
      TH1* hist = (*i)->DrawJetProfile(vars[iVar]->GetName());
      TString label = (*i)->GetPeriod(); 
      if(numData == 0) fHistManager->AddHist(hist, nLayer, label, "LP");
      else fHistManager->AddHistSame(hist, nLayer, nPad, (label += "_") +=numData, "LP");
      ++numData;
      if(replaced) (*i)->GetMultiDmFitter()->GetDstarDmFitter()
		     ->SetFitOption(dataFitOpt+=nonDataOpt);
    }

    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
//       const Int_t oldFlag = (*i)->SetFitOrNotFlag(mcFitFlag);// before 'per lumi'
      const Int_t oldFlag = (*i)->GetMultiDmFitter()->SetFitOrNotFlag(mcFitFlag);
      TH1* hist = (*i)->DrawJetProfile(vars[iVar]->GetName());
      fHistManager->AddHistSame(hist, nLayer, nPad, (*i)->GetPeriod());
      hist->SetOption("HIST");
//       (*i)->SetFitOrNotFlag(oldFlag);// before 'per lumi'
      (*i)->GetMultiDmFitter()->SetFitOrNotFlag(oldFlag);
    }
  }
  vars.Delete();
  //  fHistManager->SetNumHistsXY(4,3);
//   fHistManager->SetLogY();
  fHistManager->SetHistsMinMax(0., -1111.);
  fHistManager->Draw();
  this->Info("JetProfile", "fit flag of MC set to %d, removed '%s' from fit option of data",
	     mcFitFlag, nonDataOpt);
}


//_______________________________________________________________
 void GFDstarCompare::DrawDmShape(const char* var)
{
  // compare distributions corrected by lumi
  // NO! just normalised!
  fHistManager->Clear();
  
  Long_t numData = 0;
  TH1* histData = NULL;
  for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
    TH1* hist = (*i)->DrawDm(var);
    if(!histData) histData = hist;
    else this->Warning("DrawDmShape", "More than one data sample, normalise MC to first");
    hist->SetTitle(TString(hist->GetTitle())+= " (MC normalised to data)");
    //    GFHistManip::Normalise(hist);
    TString label = (*i)->GetPeriod(); 
//     GFDstarRunPeriods per((*i)->GetTrigger(), label);
//     hist->Scale(1./per.GetTotalLumi());
    if(numData == 0) fHistManager->AddHist(hist, 0, label);
    else fHistManager->AddHistSame(hist, 0, 0, (label += "_") +=numData);
    ++numData;
  }
  
  for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
    TH1* hist = (*i)->DrawDm(var);
//     GFHistManip::Normalise(hist);
    GFHistManip::Normalise1stTo2nd(hist, histData, "width");
//     hist->Scale(1./(*i)->GetLumi()); // is already corrected?
    fHistManager->AddHistSame(hist, 0, 0, (*i)->GetPeriod());
  }
  
  fHistManager->SetHistsLineWidth(3);
  fHistManager->Draw();
}

//_______________________________________________________________
 void GFDstarCompare::JetComposition(const char* var, Int_t flagBins)
{
  //Track,EmCl,HadCl or default(=NULL): all three
  // using fitmode 1111 in data
  // flagBins = 0: overall only, 1: overall and in pt bins, 2: only pt_bins
  TObjArray vars;
  Int_t nOverall = 0;
  if(flagBins < 2){
    if(var) {
      vars.Add(new TObjString(var));
      nOverall += 1;
    } else {
      vars.Add(new TObjString("Track"));
      vars.Add(new TObjString("EmCl"));
      vars.Add(new TObjString("HadCl"));
      nOverall += 3;
    }
  }

  Int_t nEtaBins = 0;
  Int_t nPtBins = 0;
  if(flagBins > 0){
    TH1* histJetPt = fDatas[0]->GetHist("dstar1JetPtJet", "dstarJets");
    nPtBins = (histJetPt ? histJetPt->GetNbinsY(): 0);
    for(Int_t i = 1; i <= nPtBins; ++i){
      if(var) vars.Add(new TObjString(Form("%sPt%d", var, i)));
      else {
	vars.Add(new TObjString(Form("TrackPt%d", i)));
	vars.Add(new TObjString(Form("EmClPt%d", i)));
	vars.Add(new TObjString(Form("HadClPt%d", i)));
      }
    }
    TH1* histJetEta = fDatas[0]->GetHist("dstar1JetEtaJet", "dstarJets");
    nEtaBins = (histJetEta ? histJetEta->GetNbinsY(): 0);
    for(Int_t i = 1; i <= nEtaBins; ++i){
      if(var) vars.Add(new TObjString(Form("%s%d", var, i)));
      else {
	vars.Add(new TObjString(Form("Track%d", i)));
	vars.Add(new TObjString(Form("EmCl%d", i)));
	vars.Add(new TObjString(Form("HadCl%d", i)));
      }
    }
  }

  const char* nonOpt = "M";
  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Int_t nLayer = 0;
    Int_t nPad = iVar;
    if(iVar >= nOverall){
      nLayer += (nOverall ? 1 : 0);
      nPad = iVar - nOverall;
      if(iVar >= nOverall + nPtBins*nOverall) {
	nLayer += (nPtBins ? 1 : 0);
	nPad = iVar - nOverall - nPtBins*nOverall;
      }
    }
    Long_t numData = 0;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      TString dataFitOpt = (*i)->GetMultiDmFitter()->GetDstarDmFitter()->GetFitOption();
      Bool_t repl = kFALSE;
      if(dataFitOpt.Contains(nonOpt)){
	dataFitOpt.ReplaceAll(nonOpt, NULL);
	(*i)->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(dataFitOpt);
	repl = kTRUE;
      }
      TH1* hist = (*i)->DrawJetComposition(vars[iVar]->GetName());
      hist->SetMarkerStyle(8);
      TString label = (*i)->GetPeriod(); 
      if(numData == 0) fHistManager->AddHist(hist, nLayer, label, "LP");
      else fHistManager->AddHistSame(hist, nLayer, nPad, (label += "_") +=numData, "LP");
      ++numData;
      if(repl) (*i)->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(dataFitOpt+=nonOpt);
    }

    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
// before 'per lumi'
// TString fitOp=(*i)->GetSingleMc()->GetMultiDmFitter()->GetDstarDmFitter()->GetFitOption();
      TString fitOp=(*i)->GetMultiDmFitter()->GetDstarDmFitter()->GetFitOption();
      Bool_t repl = kFALSE;
      if(fitOp.Contains(nonOpt)){
	fitOp.ReplaceAll(nonOpt, NULL);
// before 'per lumi'
// 	(*i)->GetSingleMc()->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(fitOp);
	(*i)->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(fitOp);
	repl = kTRUE;
      }
      TH1* hist = (*i)->DrawJetComposition(vars[iVar]->GetName());
      fHistManager->AddHistSame(hist, nLayer, nPad, (*i)->GetPeriod());
//       if(repl) (*i)->GetSingleMc()->GetMultiDmFitter()->GetDstarDmFitter()//before 'per lumi'
// 		 ->SetFitOption(fitOp+=nonOpt);
      if(repl) (*i)->GetMultiDmFitter()->GetDstarDmFitter()->SetFitOption(fitOp+=nonOpt);
    }
  }
  vars.Delete();

  fHistManager->SetHistsMinMax(0., kTRUE); // set minumum
  fHistManager->Draw();

  this->Info("JetComposition", "replaced fit option %s", nonOpt);
}

//_______________________________________________________________
 void GFDstarCompare::JetCompSRMean()
{
  fHistManager->Clear();
  Double_t x1, x2, y1, y2;
  fHistManager->GetLegendX1Y1X2Y2(x1, y1, x2, y2);
  cout << x1 << " " <<  x2<< " " << y1<< " " <<  y2 << endl;
  Long_t numD = 0;
  for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
    TH1* hTrackPt = (*i)->DrawJetCompSRMean("Track", "Pt");
    TH1* hTrackEta = (*i)->DrawJetCompSRMean("Track", "Eta");
    hTrackPt->SetYTitle("#bar{p_{t}(#font[32]{type})/p_{t}(jet)}");
    hTrackEta->SetYTitle(hTrackPt->GetYaxis()->GetTitle());
    hTrackPt->SetLineStyle(numD+1);
    hTrackEta->SetLineStyle(numD+1);
    TH1* hEmPt = (*i)->DrawJetCompSRMean("EmCl", "Pt");
    TH1* hEmEta = (*i)->DrawJetCompSRMean("EmCl", "Eta");
    hEmPt->SetLineStyle(numD+1);
    hEmEta->SetLineStyle(numD+1);
    hEmPt->SetLineColor(8);
    hEmEta->SetLineColor(8);
    TH1* hHadPt = (*i)->DrawJetCompSRMean("HadCl", "Pt");
    TH1* hHadEta = (*i)->DrawJetCompSRMean("HadCl", "Eta");
    hHadPt->SetLineStyle(numD+1);
    hHadEta->SetLineStyle(numD+1);
    hHadPt->SetLineColor(kMagenta);
    hHadEta->SetLineColor(kMagenta);
    TString label;// = (*i)->GetPeriod();
    TLegend *legEta = NULL, *legPt = NULL;
    if(numD == 0) {
      fHistManager->AddHist(hTrackPt, 0, "track"+label, "L");
      fHistManager->AddHist(hTrackEta, 1, "track"+label, "L");
      legPt  = fHistManager->AddLegend(0, 0, (*i)->GetPeriod(), kFALSE);
      legEta = fHistManager->AddLegend(1, 0, (*i)->GetPeriod(), kFALSE);
    } else {
      fHistManager->AddHistSame(hTrackPt, 0, 0);//, "track "+label, "L");
      fHistManager->AddHistSame(hTrackEta, 1, 0);//, "track "+label, "L");
      legPt  = new TLegend(x1-numD*(x2-x1)+0.01, y1, x2-numD*(x2-x1)+0.01,y2,(*i)->GetPeriod());
      legPt->AddEntry(hTrackPt, "track", "L");
      fHistManager->AddObject(legPt, 0, 0);
      legEta = new TLegend(x1-numD*(x2-x1)+0.01, y1, x2-numD*(x2-x1)+0.01,y2,(*i)->GetPeriod());
      legEta->AddEntry(hTrackEta, "track", "L");
      fHistManager->AddObject(legEta, 1, 0);
    }
    fHistManager->AddHistSame(hEmPt, 0, 0);//, "em "+label, "L");
    legPt->AddEntry(hEmPt, "em", "L");
    fHistManager->AddHistSame(hEmEta, 1, 0);//, "em "+label, "L");
    legEta->AddEntry(hEmEta, "em", "L");
    fHistManager->AddHistSame(hHadPt, 0, 0);//, "had "+label, "L");
    legPt->AddEntry(hHadPt, "had", "L");
    fHistManager->AddHistSame(hHadEta, 1, 0);//, "had "+label, "L");
    legEta->AddEntry(hHadEta, "had", "L");
    ++numD;
  }

  Int_t numMc = 0;
  for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
    const Int_t shift = numD+numMc;
    TH1* hTrackPt = (*i)->DrawJetCompSRMean("Track", "Pt");
    TH1* hTrackEta = (*i)->DrawJetCompSRMean("Track", "Eta");
    hTrackPt->SetLineStyle(shift+1);
    hTrackEta->SetLineStyle(shift+1);
    TH1* hEmPt = (*i)->DrawJetCompSRMean("EmCl", "Pt");
    TH1* hEmEta = (*i)->DrawJetCompSRMean("EmCl", "Eta");
    hEmPt->SetLineStyle(shift+1);
    hEmEta->SetLineStyle(shift+1);
    hEmPt->SetLineColor(8);
    hEmEta->SetLineColor(8);
    TH1* hHadPt = (*i)->DrawJetCompSRMean("HadCl", "Pt");
    TH1* hHadEta = (*i)->DrawJetCompSRMean("HadCl", "Eta");
    hHadPt->SetLineStyle(shift+1);
    hHadEta->SetLineStyle(shift+1);
    hHadPt->SetLineColor(kMagenta);
    hHadEta->SetLineColor(kMagenta);
    TString label = (*i)->GetPeriod(); 
    TLegend* legPt = new TLegend(x1-shift*(x2-x1)-0.01, y1, x2-shift*(x2-x1)-0.01,y2,
				 (*i)->GetPeriod());
    TLegend* legEta = new TLegend(x1-shift*(x2-x1)-0.01, y1, x2-shift*(x2-x1)-0.01,y2,
				  (*i)->GetPeriod());
//     if(legEta){
      fHistManager->AddHistSame(hTrackEta, 1, 0);
      fHistManager->AddHistSame(hEmEta, 1, 0);
      fHistManager->AddHistSame(hHadEta, 1, 0);
      legEta->AddEntry(hTrackEta, "track", "L");
      legEta->AddEntry(hEmEta, "em", "L");
      legEta->AddEntry(hHadEta, "had", "L");
//     } else {
//       fHistManager->AddHistSame(hTrackEta, 1, 0, "track "+label, "L");
//       fHistManager->AddHistSame(hEmEta, 1, 0, "em "+label, "L");
//       fHistManager->AddHistSame(hHadEta, 1, 0, "had "+label, "L");
//     }
//     if(legPt){
      fHistManager->AddHistSame(hTrackPt, 0, 0);
      fHistManager->AddHistSame(hEmPt, 0, 0);
      fHistManager->AddHistSame(hHadPt, 0, 0);
      legPt->AddEntry(hTrackPt, "track", "L");
      legPt->AddEntry(hEmPt, "em", "L");
      legPt->AddEntry(hHadPt, "had", "L");
//     } else {
//       fHistManager->AddHistSame(hTrackPt, 0, 0, "track "+label, "L");
//       fHistManager->AddHistSame(hEmPt, 0, 0, "em "+label, "L");
//       fHistManager->AddHistSame(hHadPt, 0, 0, "had "+label, "L");
//     }
      fHistManager->AddObject(legPt, 0, 0);
      fHistManager->AddObject(legEta, 1, 0);

    ++numMc;
  }
  const Bool_t oldDiffStyle = fHistManager->DrawDiffStyle(kFALSE);
  fHistManager->Draw();
  fHistManager->DrawDiffStyle(oldDiffStyle);
}

// //_______________________________________________________________
// void GFDstarCompare::JetPurityStability(const char* var)
// {
//   // n, eta, pt, phi, default: all
//   TObjArray vars;
//   if(var) vars.Add(new TObjString(var));
//   else {
//     vars.Add(new TObjString("pt"));
//     vars.Add(new TObjString("eta"));
//     vars.Add(new TObjString("phi"));
//     vars.Add(new TObjString("n"));
//   }
  
//   fHistManager->Clear();
//   for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
//     Long_t numMc = 0;
    
//     for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
//       TH1* histPur = (*i)->CreateJetPurity(vars[iVar]->GetName());
//       TH1* histStab = (*i)->CreateJetStability(vars[iVar]->GetName());
//       TH1* histVs = (*i)->DrawJetGenVsRec(vars[iVar]->GetName());
//       if(numMc == 0) {
// 	fHistManager->AddHist(histPur, 0, (*i)->GetPeriod());
// 	fHistManager->AddHist(histStab, 1, (*i)->GetPeriod());
//       } else {
// 	fHistManager->AddHistSame(histPur, 0, iVar, (*i)->GetPeriod());
// 	fHistManager->AddHistSame(histStab, 1, iVar, (*i)->GetPeriod());
//       }
//       fHistManager->AddHist(histVs, 2, (*i)->GetPeriod());
//       ++numMc;
//     }
//   }
//   fHistManager->AddHistsOption("COLZ", 2);
//   vars.Delete();

//   fHistManager->Draw();
// }

//_______________________________________________________________
 void GFDstarCompare::DsJetPurityStability(const char* var, Bool_t hadCorr)
{
  //Dphi Deta xGam PtJet PtDs EtaJet default: all
  TObjArray vars;
  if(var) vars.Add(new TObjString(var));
  else {
    vars.Add(new TObjString("Dphi"));
    vars.Add(new TObjString("Deta"));
    vars.Add(new TObjString("xGam"));
    vars.Add(new TObjString("PtJet"));
    vars.Add(new TObjString("EtaJet"));
    if(!hadCorr) {
      vars.Add(new TObjString("PtDs"));
      vars.Add(new TObjString("CosTh"));
      vars.Add(new TObjString("MDsJet"));
      vars.Add(new TObjString("PtDsJet"));
      vars.Add(new TObjString("ChFrac"));
      vars.Add(new TObjString("DaughJet"));
      vars.Add(new TObjString("PhiJet"));
      vars.Add(new TObjString("DeltaR"));
      vars.Add(new TObjString("All"));
    }
  }
  
  fHistManager->Clear();
  for(Int_t iVar = 0; iVar < vars.GetEntriesFast(); ++iVar){
    Long_t numMc = 0;
    
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* histPur = (hadCorr ? (*i)->CreateDsJetPurityHadCorr(vars[iVar]->GetName()) 
		              : (*i)->CreateDsJetPurity(vars[iVar]->GetName()));
      TH1* histStab = (hadCorr ? (*i)->CreateDsJetStabilityHadCorr(vars[iVar]->GetName()) 
		               : (*i)->CreateDsJetStability(vars[iVar]->GetName()));
      //      TH1* histVs = (*i)->DrawJetGenVsRec(vars[iVar]->GetName());
      if(numMc == 0) {
	fHistManager->AddHist(histPur, 0, (*i)->GetPeriod());
 	fHistManager->AddHist(histStab, 1, (*i)->GetPeriod());
      } else {
	fHistManager->AddHistSame(histPur, 0, iVar, (*i)->GetPeriod());
 	fHistManager->AddHistSame(histStab, 1, iVar, (*i)->GetPeriod());
      }
//       fHistManager->AddHist(histVs, 2, (*i)->GetPeriod());
      ++numMc;
    }
  }
//   fHistManager->AddHistsOption("COLZ", 2);
  vars.Delete();

  fHistManager->AddHistsOption("E");
  fHistManager->Draw();
}


// //_______________________________________________________________
// TH1* GFDstarCompare::JetResEta(Int_t etaBin)
// {
//   // jet resolution and rec vs gen in eta bin (etaBin == 0) -> all bins
//   // n, eta, pt, phi, default: all

//   Int_t firstBin = etaBin, lastBin = etaBin;
// //   TH1 *hMean1 = NULL, *hMean2 = NULL, *hSigma1 = NULL, *hSigma2 = NULL;

//   if(etaBin == 0){
//     firstBin = 1;
//     TH1* h = fMonteCarlo[0]->CreateJetPurity("eta");
//     lastBin = h->GetNbinsX();
// //     const Double_t xMin = h->GetXaxis()->GetXmin();
// //     const Double_t xMax =  h->GetXaxis()->GetXmax();
// //     hMean1 = new TH1F("mean1", "#mu_{1} of di-gaus", lastBin, xMin, xMax);
// //     hMean2 = new TH1F("mean2", "#mu_{2} of di-gaus", lastBin, xMin, xMax);
// //     hSigma1 = new TH1F("sigma1", "#sigma_{1} of di-gaus", lastBin, xMin, xMax);
// //     hSigma2 = new TH1F("sigma2", "#sigma_{2} of di-gaus", lastBin, xMin, xMax);
//   }

//   fHistManager->Clear();

// //   TF1* func = GFFuncs::CreateDiGaus1Mean("digaus1m");
// //   TF1* func = GFFuncs::CreateDiGaus("digaus1m", -.5, .5);
// //   TF1* func = GFFuncs::CreateGaus("mygaus", -.5, .5);

//   const Bool_t addLegend = (fMonteCarlo.size() > 1);


//   for(Int_t bin = firstBin; bin <= lastBin; bin++){

//     Long_t numMc = 0;
//     for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
//       TH1* histResRel = (*i)->DrawJetRes("pt", "Rel", bin);
// //       TH1* histRes = (*i)->DrawJetRes("pt", "", bin);
//       if(numMc == 0) {
//  	fHistManager->AddHist(histResRel, 0, addLegend ? (*i)->GetPeriod() : NULL);
// // 	fHistManager->AddHist(histResRel, 0);
// 	// 	histResRel->Fit(func, "WIR");
// // 	hMean1->SetBinContent(bin, func->GetParameter(1));
// // // 	hMean1->SetBinError(bin, func->GetParError(1));
// // 	hMean2->SetBinContent(bin, func->GetParameter(4));
// // // 	hMean2->SetBinError(bin, func->GetParError(4));
// // 	hSigma1->SetBinContent(bin, func->GetParameter(2));
// // // 	hSigma1->SetBinError(bin, func->GetParError(2));
// // 	hSigma2->SetBinContent(bin, func->GetParameter(5));
// // // 	hSigma2->SetBinError(bin, func->GetParError(5));
// //   	fHistManager->AddHist(histRes, bin, (*i)->GetPeriod());
//       } else {
// // 	continue;
// 	fHistManager->AddHistSame(histResRel, 0, bin-firstBin, addLegend ? (*i)->GetPeriod() : NULL);
// //   	fHistManager->AddHistSame(histRes, bin, 1, addLegend ? (*i)->GetPeriod() : NULL);
//       }
// //       TH1* histVs = (*i)->DrawJetGenVsRec("pt", bin);
// //       fHistManager->AddHist(histVs, lastBin+1, (*i)->GetPeriod());
//       ++numMc;
//     }
//   }

// //   delete func;

// //   if(hMean1)  fHistManager->AddHist(hMean1, 1);
// //   if(hMean2)  fHistManager->AddHist(hMean2, 1);
// //   if(hSigma1) fHistManager->AddHist(hSigma1,1);
// //   if(hSigma2) fHistManager->AddHist(hSigma2,1);
// //   fHistManager->AddHistsOption("COLZ", lastBin+1);

//   fHistManager->Draw();
//   TCanvas *can = fHistManager->GetCanvas(0,0);
//   for(Int_t i = 1; i <= 8; ++i){
//     can->cd(i);
//     TLine* l = new TLine(0,0,0,1500); 
//     l->Draw();
//   }
//   return NULL;
// }

//_______________________________________________________________
 void GFDstarCompare::Ds1Jet(const char* forwBack)
{
  // "Forw", "Back", ""
  fHistManager->Clear();
  
  TObjArray var;
  var.Add(new TObjString("EtaDs"));
  var.Add(new TObjString("PtDs"));
  var.Add(new TObjString("EtaJet"));
  var.Add(new TObjString("PtJet"));
  var.Add(new TObjString("Deta"));
  var.Add(new TObjString("Dphi"));
  var.Add(new TObjString("Dpt"));

  for(Int_t iVar = 0; iVar < var.GetEntriesFast(); ++iVar){
    Long_t numData = 0;
    TH1* histData = NULL;
    for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
      TH1* hist = (*i)->DrawDmDs1Jet(var[iVar]->GetName(), forwBack);
      if(!histData) histData = hist;
      else this->Warning("Ds1Jet", "More than one data sample, normalise MC to first");
      hist->SetTitle(TString(hist->GetTitle())+= " (MC normalised to data)");
      //    GFHistManip::Normalise(hist);
      TString label = (*i)->GetPeriod(); 
      //     GFDstarRunPeriods per((*i)->GetTrigger(), label);
      //     hist->Scale(1./per.GetTotalLumi());
      if(numData == 0) fHistManager->AddHist(hist, 0, label);
      else fHistManager->AddHistSame(hist, 0, iVar, (label += "_") +=numData);
      ++numData;
    }
  
    for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
      TH1* hist = (*i)->DrawDmDs1Jet(var[iVar]->GetName(), forwBack);
      //     GFHistManip::Normalise(hist);
      GFHistManip::Normalise1stTo2nd(hist, histData, "width");
      //     hist->Scale(1./(*i)->GetLumi()); // is already corrected?
      fHistManager->AddHistSame(hist, 0, iVar, (*i)->GetPeriod());
    }
  }
  
//   fHistManager->SetHistsLineWidth(3);
  fHistManager->Draw();
  var.Delete();
}

//_______________________________________________________________
 void GFDstarCompare::TrackErr(Option_t *option)
{
// option may contain "rel" [default] to see sigma(p)/p instead of sigma(p)

  fHistManager->Clear();

  const TString opt(option);
  TString name("sigmaP");
  if (opt.Contains("rel", TString::kIgnoreCase)) name += "rel";

  TH1 *profPtThetaData = NULL;
  for (DataIter i = fDatas.begin(); i != fDatas.end(); ++i) {
    TH1 *h = (*i)->GetHist(name, "trackErr");
    Double_t scale = h->Integral();
    if (!scale) {
      scale = 1.;
      this->Warning("TrackErr", "empty hist data");
    } else scale = 1./scale;
    h->Scale(scale);
    fHistManager->AddHistSame(h, 0, 0, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsPt", "trackErr");
    h->Scale(scale);
    fHistManager->AddHistSame(h, 1, 0, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsTheta", "trackErr");
    h->Scale(scale);
    fHistManager->AddHistSame(h, 1, 1, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsPtTheta", "trackErr");
//     h->Scale(scale);
    if (!profPtThetaData) profPtThetaData = h;
    fHistManager->AddHist(h, 2, (*i)->GetPeriod());
  }
  for (McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i) {
    TH1 *h = (*i)->GetHist(name, "trackErr");
    Double_t scale = h->Integral();
    if (!scale) {
      scale = 1.;
      this->Warning("TrackErr", "empty hist mc");
    } else scale = 1./scale;
    h->Scale(scale);
    fHistManager->AddHistSame(h, 0, 0, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsPt", "trackErr");
    h->Scale(scale);
    fHistManager->AddHistSame(h, 1, 0, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsTheta", "trackErr");
    h->Scale(scale);
    fHistManager->AddHistSame(h, 1, 1, (*i)->GetPeriod());
    h = (*i)->GetHist(name + "vsPtTheta", "trackErr");
//     h->Scale(scale);
    fHistManager->AddHist(h, 2, (*i)->GetPeriod());
    if (profPtThetaData) {
      TH1 *rel = 
	GFHistManip::CreateAnalog(profPtThetaData, 
				  Form("%s%s", h->GetName(), (*i)->GetPeriod()),
				  kFALSE, kTRUE); // no err, but fill
      rel->Divide(profPtThetaData, h);
      rel->SetTitle(Form("%s, data divided by %s",
			 rel->GetTitle(), (*i)->GetPeriod()));
      fHistManager->AddHist(rel, 3);
      TH1 *hRat = new TH1F(Form("ratio%s", (*i)->GetPeriod()), rel->GetTitle(),
			   40,0.5,1.5);
      hRat->SetXTitle(h->GetTitle());
      for (Int_t iX = 1; iX <= rel->GetNbinsX(); ++iX) {
	for (Int_t iY = 1; iY <= rel->GetNbinsY(); ++iY) {
	  const Float_t ratio = rel->GetBinContent(iX, iY);
	  if (ratio) hRat->Fill(ratio);
	}
      }
      fHistManager->AddHistSame(hRat, 4, 0, (*i)->GetPeriod());
    }
  }
  
  fHistManager->SetLogY(0, kTRUE);
  fHistManager->AddHistsOption("BOX", 1);
  fHistManager->AddHistsOption("surf1", 2);
  fHistManager->SetHistsMinMax(.1, kFALSE, 2);
  if (fHistManager->GetNumLayers() > 4) {
    fHistManager->AddHistsOption("surf3z", 3);
    fHistManager->SetHistsMinMax(2., kFALSE, 3);
    fHistManager->SetHistsMinMax(.6, kTRUE, 3);
    fHistManager->AddHistsOption("histsames", 4);
    fHistManager->GetHistsOf(4, 0)->First()->SetOption("");
    fHistManager->SetLogY(4, kTRUE);
  }
  fHistManager->Draw();
  this->Info("TrackErr", "Pythia hist might be wrongly merged!");
}


//_______________________________________________________________
 void GFDstarCompare::SetTrigger(Int_t subTr)
{
  for(DataIter i = fDatas.begin(); i != fDatas.end(); ++i){
    (*i)->SetTrigger(subTr);
  }

  for(McIter i = fMonteCarlo.begin(); i != fMonteCarlo.end(); ++i){
    (*i)->SetTrigger(subTr);
  }
}

//_______________________________________________________________
 GFDstarHistsAnalysisData* GFDstarCompare::GetHistsAnalysisData(Int_t i)
{
  return fDatas[i];
}

//_______________________________________________________________
 GFDstarHistsAnalysisMc* GFDstarCompare::GetHistsAnalysisMc(Int_t i)
{
  return fMonteCarlo[i];
}


ROOT page - Class index - Top of the page

This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.