// Author: Gero Flucke <mailto:flucke@mail.desy.de>
//____________________________________
// GFDstarXSections
//   Author:      Gero Flucke
//   Date:        May 14th, 2002
//   last update: $Date: 2005/11/16 00:39:02 $  
//   by:          $Author: flucke $
//

#include <iostream> // debugging only
using namespace std;
#include <limits.h>

#include <TROOT.h>
#include <TFile.h>
#include <TPad.h>
#include <TCanvas.h>
#include <TH1.h>
#include <TGraphAsymmErrors.h>
#include <TLegend.h>
#include <TLegendEntry.h>
// #include <TH2.h>
#include <TString.h>
#include <TObjString.h>
// #include <TAxis.h>
// #include <TObjArray.h>

#include "GFAnalysis/GFDstarXSections.h"
#include "GFAnalysis/GFDstarHistsAnalysisData.h"
#include "GFAnalysis/GFDstarHistsAnalysisMc.h"
#include "GFAnalysis/GFAxisVariables.h"
#include "GFTools/GFFitHistFrac.h"

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

ClassImp(GFDstarXSections)

 GFDstarXSections::GFDstarXSections(const char* filePrefix, const char* period)
  : fTrigger(83), fInclXsec(2), fMcFactors(0)
{

  fAnalysisData = new GFDstarHistsAnalysisData(filePrefix, (period ? period : "All"));
  fAnalysisMc   = new GFDstarHistsAnalysisMc(filePrefix, 
					     (period ? Form("Pythia%s", period) : "Pythia"));


  fAnalysisOtherMc = new GFDstarHistsAnalysisMc(filePrefix, 
						(period ? Form("Cascade%s", period):"Cascade"));

  this->Init();
}

//____________________________________
 GFDstarXSections::GFDstarXSections(GFDstarHistsAnalysisData* data, 
				   GFDstarHistsAnalysisMc* m1,
				   GFDstarHistsAnalysisMc* m2) :
  fTrigger(83), fInclXsec(2), fMcFactors(0)
{
  // pointers are adopted!
  fAnalysisData = data;
  fAnalysisMc = m1;
  fAnalysisOtherMc = m2;

  if(!fAnalysisMc || !fAnalysisData){
    this->Error("GFDstarXSections", "data or mc missing!");
  }
  fTrigger = fAnalysisData->GetTrigger();
  if(fTrigger != fAnalysisMc->GetTrigger()){
    this->Warning("GFDstarXSections", "trigger %d in data, MC was %d, use %d", fTrigger, 
		  fAnalysisMc->GetTrigger(), fTrigger);
    fAnalysisMc->SetTrigger(fTrigger);
  }
  if(fAnalysisOtherMc && fTrigger != fAnalysisOtherMc->GetTrigger()){
        this->Warning("GFDstarXSections", "trigger %d in data, other MC was %d, use %d", 
		      fTrigger, fAnalysisOtherMc->GetTrigger(), fTrigger);
	fAnalysisOtherMc->SetTrigger(fTrigger);
  }

  this->Init();
}

//____________________________________
 void GFDstarXSections::Init()
{
  fHistManager  = new GFHistManager;

  this->SetLumi(); // after fAnalysisData is set!
  fAcceptFromMc.Set(2);
  fAcceptFromMc[1] = kTRUE; //S83
  fAcceptFromMc[0] = kTRUE;  //S84
  fEffTrigFromMc = fAcceptFromMc;

  fReflectCorrection = 1. - 0.035;// from S. Hengstmann's thesis
  fEffDedx           = 0.98;

  fAnalysisData->SetBatch();
  fAnalysisMc->SetBatch();
  Int_t nMc = (fAnalysisMc ? 1 : 0);
  if(fAnalysisOtherMc) {
    fAnalysisOtherMc->SetBatch();
    ++nMc;
  }
  fMcFactors.Set(nMc);

  this->SetHistPointerNull();
}

//____________________________________
 GFDstarXSections::~GFDstarXSections()
{
  if(fHistManager) delete fHistManager;

  if(fAnalysisData) delete fAnalysisData;
  if(fAnalysisMc) delete fAnalysisMc;
  if(fAnalysisOtherMc) delete fAnalysisOtherMc;

  if(fHistNumDstar) delete fHistNumDstar;
  if(fHistTrigEff) delete fHistTrigEff;
  if(fHistRecEff)  delete fHistRecEff;
  if(fHistAccept)  delete fHistAccept;
  if(fHistData1) delete fHistData1;
  if(fHistData2) delete fHistData2;
  if(fHistData3) delete fHistData3;
  if(fHistMc1) delete fHistMc1;
  if(fHistMc2) delete fHistMc2;
  if(fHistMc3) delete fHistMc3;
}

//____________________________________
 TH1* GFDstarXSections::XSectionNorm(const char *var)
{
  const Bool_t wasBatch = fHistManager->SetBatch(kTRUE);
  TH1 *hist = this->XSection(var);
  if (hist) {
    hist = this->XSectionNormInt(hist, "D*", var);
  }
  fHistManager->SetBatch(wasBatch);
  fHistManager->Draw();
  return hist;
}

//____________________________________
 TH1* GFDstarXSections::XSectionNorm(const char *var1, const char *var2,
				    Int_t firstBin2, Int_t lastBin2)
{
  const Bool_t wasBatch = fHistManager->SetBatch(kTRUE);
  TH1 *hist = this->XSection(var1, var2, firstBin2, lastBin2);
  if (hist) {
    hist = this->XSectionNormInt(hist, "D*", var1);
  }
  fHistManager->SetBatch(wasBatch);
  fHistManager->Draw();
  return hist;
}

//____________________________________
 TH1* GFDstarXSections::XSectionDsJetNorm(const char *var, const char* forwBack)
{
  const Bool_t wasBatch = fHistManager->SetBatch(kTRUE);
  TH1 *hist = this->XSectionDsJet(var, forwBack);
  if (hist) {
    const TString sample(Form("D* + oth. jet%s", forwBack));
    hist = this->XSectionNormInt(hist, sample.Data(), var);
  }
  fHistManager->SetBatch(wasBatch);
  fHistManager->Draw();
  return hist;
}

//____________________________________
 TH1* GFDstarXSections::XSectionDsDiJetNorm(const char* var)
{
  const Bool_t wasBatch = fHistManager->SetBatch(kTRUE);
  TH1 *hist = this->XSectionDsDiJet(var);
  if (hist) {
    hist = this->XSectionNormInt(hist, "D* + 2 jet", var);
  }
  fHistManager->SetBatch(wasBatch);
  fHistManager->Draw();
  return hist;
}

//____________________________________
 TH1* GFDstarXSections::XSectionNormInt(TH1 *hist, const char *sample, const char *var)
{
// FIXME: MC?
  if (!hist) return NULL;
  TH1 *hNew = static_cast<TH1*>(hist->Clone(Form("%sNorm", hist->GetName())));
  GFHistManip::ReCorrectForBinWidth(hNew);
  Double_t totXSec = 0., totXSecErr2 = 0.;// sum of bins, including over-/underflow
  const Int_t overflowBin = hNew->GetNbinsX() + 1;
  for (Int_t iBin = 0; iBin <= overflowBin; ++iBin) {
    totXSec += hNew->GetBinContent(iBin);
    totXSecErr2 += hNew->GetBinError(iBin) * hNew->GetBinError(iBin);
  }

  TArrayD tmp(2);
  tmp[0] = totXSec;
  tmp[1] = TMath::Sqrt(totXSecErr2);
  fAnalysisData->PrintPlusMinus(tmp, Form("total X-section as sum of previous"));

  const Double_t totXSec2 = totXSec*totXSec;

  for (Int_t iBin = 0; iBin <= overflowBin; ++iBin) { // FIXME: skip over-/underflow?
    const Double_t cont = hNew->GetBinContent(iBin);
    const Double_t err2 = hNew->GetBinError(iBin) * hNew->GetBinError(iBin);
    hNew->SetBinContent(iBin, cont/totXSec);
// cf. log book November 15th, 2005:
    const Double_t newErr2 = 1./totXSec2 * (err2 + cont/totXSec *
					    (totXSecErr2 * cont/totXSec -2. * err2));
    hNew->SetBinError(iBin, TMath::Sqrt(newErr2));
  }
  GFHistManip::CorrectForBinWidth(hNew);
  fHistManager->AddHist(hNew, fHistManager->GetNumLayers());

  TString yTit(Form("1/#sigma_{%s} %s", sample, hNew->GetYaxis()->GetTitle()));
  TString tit(Form("1/#sigma_{%s} %s", sample, hNew->GetTitle()));
  const char *unit = this->GetUnit(GFAxisVariables::GetAxisLabel(var));
  if (unit) {
    yTit.ReplaceAll("nb", "1");
    tit.ReplaceAll(Form("[%s]", unit), NULL);
  } else {
    yTit.ReplaceAll("[nb]", NULL);
    tit.ReplaceAll("[nb]", NULL);
  }
  hNew->SetTitle(Form("%s;%s;%s", tit.Data(), hNew->GetXaxis()->GetTitle(), yTit.Data()));

  GFHistManip::PrintTable(hNew, "previous now normalised"); 
  return hNew;
}


//____________________________________
 TH1* GFDstarXSections::XSectionDsJet(const char* var, const char* forwBack)
{
  fHistNumDstar = fAnalysisData->DrawDmDs1Jet(var, forwBack);
  fHistTrigEff  = fAnalysisMc->CreateTrigEffDsJet(var, forwBack);//FIMXE ?: CreateTrigEff(var);
  fHistRecEff   = fAnalysisMc->CreateRecEffDs1Jet(var, forwBack);
  fHistAccept   = fAnalysisMc->CreateAcceptDsJet(var, forwBack);// FIXME ?: always from MC?
  fHistMc1      = fAnalysisMc->XSectionDsJet(var, forwBack);
  fHistMc2      = fAnalysisOtherMc ? fAnalysisOtherMc->XSectionDsJet(var, forwBack) : NULL;
  
  TString pro(Form("D*%s + 1 jet", forwBack));
  return this->XSectionInt(GFAxisVariables::GetAxisLabel(var), pro.Data());
}

//____________________________________
 TH1* GFDstarXSections::XSectionDsDiJet(const char* var)
{
  fHistNumDstar = fAnalysisData->DrawDmDsDiJet(var);
  fHistTrigEff  = fAnalysisMc->CreateTrigEffDsDiJet(var);//FIMXE ?: CreateTrigEff(var);
  fHistRecEff   = fAnalysisMc->CreateRecEffDsDiJet(var);
  fHistAccept   = fAnalysisMc->CreateAcceptDsDiJet(var);// FIXME ?: always from MC?
  fHistMc1      = fAnalysisMc->XSectionDsDiJet(var);
  fHistMc2      = fAnalysisOtherMc ? fAnalysisOtherMc->XSectionDsDiJet(var) : NULL;
  
  return this->XSectionInt(GFAxisVariables::GetAxisLabel(var), "D* dijet");
}

//____________________________________
 TH1* GFDstarXSections::XSection(const char* var1, const char* var2, Int_t firstBin2,
				Int_t lastBin2)
{
  fHistNumDstar = fAnalysisData->DrawDm(var1, var2, firstBin2, lastBin2, kTRUE);
  fHistTrigEff  = fAnalysisMc->CreateTrigEff(var1, var2, firstBin2, lastBin2);
  fHistRecEff   = fAnalysisMc->CreateRecEff(var1, var2, firstBin2, lastBin2);
  fHistAccept   = fAnalysisMc->CreateAccept(var1, var2, firstBin2, lastBin2);

  fHistMc1      = fAnalysisMc->XSection(var1, var2, firstBin2, lastBin2);
  fHistMc2      = (fAnalysisOtherMc ? 
		   fAnalysisOtherMc->XSection(var1, var2, firstBin2, lastBin2) : NULL);
  
  TH1* tmp = fHistMc1;
  TH1* hXSec = this->XSectionInt(GFAxisVariables::GetAxisLabel(var1), "D*");
  if(hXSec) hXSec->SetTitle(tmp->GetTitle());
  return hXSec;
}

//____________________________________
 TH1* GFDstarXSections::XSection(const char* var, Bool_t acceptOnMc)
{
  fHistNumDstar = fAnalysisData->DrawDm(var);
  fHistTrigEff  = this->CreateTrigEff(var);
  fHistRecEff   = fAnalysisMc->CreateRecEff(var);
  if(acceptOnMc) fHistAccept = NULL;
  else{
    fHistAccept   = this->AcceptFromMc() ? fAnalysisMc->CreateAccept(var)
      : fAnalysisData->DrawAccept(var);// FIXME ???
  }

  fHistMc1      = fAnalysisMc->XSection(var, acceptOnMc);
  fHistMc2      = fAnalysisOtherMc ? fAnalysisOtherMc->XSection(var, acceptOnMc) : NULL;

  return this->XSectionInt(GFAxisVariables::GetAxisLabel(var), "D*");
}

//____________________________________
 TH1* GFDstarXSections::XSectionWgammaP(Bool_t acceptOnMc)
{
  this->InitXSecWgammaP(acceptOnMc);
  return this->XSectionInt("W_{#gammap}", "D*");
}

//____________________________________
 void GFDstarXSections::InitXSecWgammaP(Bool_t acceptOnMc)
{
  // set internal hist pointer for XSectioWgammaP...()
  fHistNumDstar = fAnalysisData->DrawDm("wGammaP");//WgammaP();
  fHistTrigEff  = this->CreateTrigEff("wGammaP");
  fHistRecEff   = fAnalysisMc->CreateRecEff("wGammaP");
  if(acceptOnMc) fHistAccept = NULL;
  else {
    fHistAccept   = this->AcceptFromMc() ? fAnalysisMc->CreateAccept("wGammaP")
                                         : fAnalysisData->DrawAccept("wGammaP");//old: WgammaP
  }
  // Has to be the last for XSectionWgammaPBoth:
  fHistMc1 = fAnalysisMc->XSectionWgammaP(acceptOnMc);
  fHistMc2 = fAnalysisOtherMc ? fAnalysisOtherMc->XSection("wGammaP", acceptOnMc) : NULL;
}


//____________________________________
 TH1* GFDstarXSections::XSectionInt(const char* variable, const char* process, 
				   Bool_t drawWrite)
{
  //  if(drawWrite) fHistManager->Clear();
  fHistManager->Clear(); //FIXME

  if(!fHistNumDstar || !fHistRecEff || !fHistTrigEff || !fHistMc1) {
    return NULL;// fHistAccept or fHistMc2 may be missing
  }

  Int_t numBins = fHistNumDstar->GetNbinsX();
  TArrayD relErrorNumDstar(numBins);
  for(Int_t i = 0; i < numBins; ++i){
    relErrorNumDstar[i] = fHistNumDstar->GetBinError(i+1);
    if(fHistNumDstar->GetBinContent(i+1)) {
      relErrorNumDstar[i] /= fHistNumDstar->GetBinContent(i+1);
    }
  }

  GFHistManip::Scale(fHistNumDstar, fReflectCorrection);

  TString name("xSection");
  fAnalysisData->MakeStringName(name += variable);

  TH1* temp = static_cast<TH1*>(fHistTrigEff->Clone(name+"Temp"));
  temp->Multiply(fHistTrigEff, fHistRecEff, 
		 fAnalysisMc->GetBranchRatio(), this->GetLumi());

  if(!fHistAccept){
    this->Warning("XSection", "using mode ignoring acceptance!");
  } else if(fTrigger == 84){
    if(this->AcceptFromMc()){
      temp->Multiply(temp, fHistAccept);
    } else {
      temp->Scale(fAnalysisData->TotalAccept().At(0));
      this->Warning("XSection", 
		    "ETAG44: ignore acceptance histogram, use overall factor from data");
    }
//   } else if(fTrigger == 83){
//     TString varStr(variable);
//     if(varStr == "W_{#gammap}" || (varStr.Contains("eta") && !varStr.Contains("-"))){
//       // w_gp, eta, but not delta eta...
//       this->Info("XSection", "napplying acceptance histogram!n");
//       temp->Multiply(temp, fHistAccept);
//     } else {
//       this->Warning("XSection", "using overall factor for acceptance"); 
//       Double_t totalAccept = this->AcceptFromMc() ? fAnalysisMc->TotalAccept().At(0)
//                                                   : fAnalysisData->TotalAccept().At(0);
//       temp->Scale(totalAccept);
//     }
  } else {
    temp->Multiply(temp, fHistAccept);
  }

  TH1* histXSection = static_cast<TH1*>(fHistNumDstar->Clone(name+=fTrigger));
  const Double_t l4Eff = fAnalysisData->L4Efficiency().At(0);
  histXSection->Divide(histXSection, temp, 1., l4Eff*fEffDedx);
  delete temp;

  for(Int_t i = 0; i < numBins; ++i){
    histXSection->SetBinError(i+1, relErrorNumDstar[i] * histXSection->GetBinContent(i+1));
  }
  TString title(Form("d#sigma_{vis} (ep#rightarrow %s X) / d", process));
  title+=variable;//)+=" [pb]";
  histXSection->SetTitle(title);

  TString titleShort(Form("d#sigma / d%s [nb", variable));
  const char *unit = this->GetUnit(titleShort);//NULL;
//   if(titleShort.Contains("[GeV]")) unity = "GeV";
//   else if(titleShort.Contains("[#circ]")) unity = "#circ";
  if(unit){
    titleShort.ReplaceAll(Form("[%s]", unit), NULL); // pb
    (titleShort += "/") += unit;
  }
  titleShort += "]";
  histXSection->SetYTitle(titleShort);
  histXSection->SetMarkerStyle(8);//kFullDotLarge);

  fHistManager->AddHist(fHistNumDstar, 0);
  fHistManager->AddHist(fHistTrigEff, 0);
  fHistManager->AddHist(fHistRecEff, 0);
  fHistManager->AddHist(fHistAccept, 0); // if NULL pointer: ignored!
  TString per;// = "H1 ";
  per += fAnalysisData->GetPeriod();
  fHistManager->AddHist(histXSection, 1, per.Data(), "lp");
  if(!TString(variable).Contains("p_{t}")) histXSection->SetMinimum(0.); // mmhhh: fixme?
  else {
    histXSection->SetMinimum(-1111.);
    fHistManager->SetLogY(1);
  }
  
  fMcFactors[0] = 1.;//drawWrite ? this->NormaliseToData(fHistMc1, fAnalysisMc) : 1.;
  //  fHistMc1->SetOption("HIST");
  fHistMc1->SetMaximum(-1111.);//fHistMc1->GetMaximum() * factor);
//   fHistManager->AddHistSame(fHistMc1, 1, 0, TString("PYTHIA #times ")+=factor);
  fHistManager->AddHistSame(fHistMc1, 1, 0, Form("%s #times %f", 
						 fAnalysisMc->GetPeriod(), fMcFactors[0]));
  fHistManager->AddLegend(1, 0, fAnalysisData->GetSampleName());

  if(fHistMc2 && fAnalysisOtherMc){
    fMcFactors[1] = 1.;//this->NormaliseToData(fHistMc2, fAnalysisOtherMc);
//     fHistMc2->SetOption("HIST");
    fHistMc2->SetMaximum(-1111.);//fHistMc2->GetMaximum() * factor);
//     fHistManager->AddHistSame(fHistMc2, 1, 0, TString("Cascade #times ")+=factor2);
    fHistManager->AddHistSame(fHistMc2, 1, 0, 
			      Form("%s #times %f",fAnalysisOtherMc->GetPeriod(),fMcFactors[1]));
  }

//   GFHistArray* singleMcs = fAnalysisMc->GetHistsTmp();
//   for(Int_t i = 0; i < singleMcs->GetEntriesFast(); ++i){
//     (*singleMcs)[i]->Scale(fMcFactors[1]);
//     (*singleMcs)[i]->SetMaximum(-1111.);//(*singleMcs)[i]->GetMaximum() * factor);
//     fHistManager->AddHistSame((*singleMcs)[i], 1, 0, fAnalysisMc->GetSingleMc(i)->GetPeriod());
//   }


//   GFHistArray* mcFitToData = this->FitFractions(histXSection, fAnalysisMc->GetHistsTmp());
//   fHistManager->AddHistsSame(mcFitToData, 1, 1);
//   delete mcFitToData;  mcFitToData = NULL;

  if(drawWrite){
    fHistManager->Draw();
    this->SetHistPointerNull();
  }
  title = Form("%s X-section in bins of ", process);
  fAnalysisData->PrintTable(histXSection, title += variable); 
  return histXSection;
}


//____________________________________
 TArrayD GFDstarXSections::XSectionTot(Bool_t ignoreAccept) const
{
  if(ignoreAccept || !fInclXsec[0]){

    const TArrayD numDstar = fAnalysisData->DrawDm();
    const TArrayD trigEff = this->TotalTrigEff();
    const TArrayD recEff  = fAnalysisMc->TotalRecEff();
    const TArrayD l4Eff   = fAnalysisData->L4Efficiency();
    Double_t accept = 1.;
    if(!ignoreAccept) {
      accept = this->AcceptFromMc() ? fAnalysisMc->TotalAccept().At(0) 
	: fAnalysisData->TotalAccept().At(0);
    }
    
    fInclXsec = this->XSectionTot(numDstar, trigEff[0], recEff[0], accept, l4Eff[0], "D*");
  }

  return fInclXsec;
}

 TArrayD GFDstarXSections::XSectionTot(const TArrayD& numDstar, Double_t trigEff, 
				      Double_t recEff, Double_t accept, Double_t l4Eff,
				      const char* what) const
{
  TArrayD xSectAndError(2);
  const Double_t nominator = fAnalysisMc->GetBranchRatio() * trigEff * recEff
    * this->GetLumi() * accept * l4Eff * fEffDedx;
  xSectAndError[0] = nominator ? numDstar[0]/nominator : 0.;
  xSectAndError[1] = numDstar[0] ? xSectAndError[0] * numDstar[1]/numDstar[0] : 0.;

  // reflection correction
  xSectAndError[0] *= fReflectCorrection;
  xSectAndError[1] *= fReflectCorrection;
  
  fAnalysisData->PrintPlusMinus(xSectAndError, Form("%s: Total X-section", what));

  return xSectAndError;
}

//____________________________________
 TArrayD GFDstarXSections::XSectionDsJetTot(const char* forwBack) const
{
  const TArrayD numDstar = fAnalysisData->DrawDmDsJet(forwBack);
  const TArrayD trigEff = fAnalysisMc->TotalTrigEffDsJet(forwBack);
  const TArrayD recEff  = fAnalysisMc->TotalRecEffDsJet(forwBack);
  const TArrayD l4Eff   = fAnalysisData->L4Efficiency();
  const Double_t accept = fAnalysisMc->TotalAcceptDsJet(forwBack).At(0);
    
  return this->XSectionTot(numDstar, trigEff[0], recEff[0], accept, l4Eff[0], 
			   Form("D* %s + 1 jet", forwBack));
}

//____________________________________
 TArrayD GFDstarXSections::XSectionDsDiJetTot() const
{
  const TArrayD numDstar = fAnalysisData->DrawDmDsDiJet();
  const TArrayD trigEff = fAnalysisMc->TotalTrigEffDsDiJet();
  const TArrayD recEff  = fAnalysisMc->TotalRecEffDsDiJet();
  const TArrayD l4Eff   = fAnalysisData->L4Efficiency();
  const Double_t accept = fAnalysisMc->TotalAcceptDsDiJet().At(0);
    
  return this->XSectionTot(numDstar, trigEff[0], recEff[0], accept, l4Eff[0], "D* dijet");
}


//____________________________________
 TH1* GFDstarXSections::XSectionTotHist(Bool_t draw, Bool_t ignoreAccept)
{
  if(draw) fHistManager->Clear();

  TArrayD xSectAndError = this->XSectionTot(ignoreAccept);
  TString name("sigmaTot");
  TString title("#sigma_{vis} ep #rightarrow D* + X [nb]"); // pb


  TH1* hTemp = fAnalysisData->GetHist2D("wGammaP");
  TH1* histX = new TH1F(name+=fTrigger, 
			title+=fAnalysisData->GetSampleName(), 
			1, hTemp->GetYaxis()->GetXmin(), hTemp->GetYaxis()->GetXmax());
  histX->SetBinContent(1, xSectAndError[0]);
  histX->SetBinError(1, xSectAndError[1]);
  
  if(draw){
    fHistManager->AddHist(histX);
    fHistManager->Draw();
  }
  return histX;
}


//____________________________________
 TH1* GFDstarXSections::CreateTrigEff(const char* var)
{
  TH1* hTr19 = NULL;
  if(this->GetTrigger() == 83){
    hTr19 = this->EffTrigFromMc() ? 
      fAnalysisMc->CreateTrigEff(var) : fAnalysisData->CreateTrigEff(var);
  } else if(this->GetTrigger() == 84){
    this->Warning("CreateTrigEff", "Mixing MC per bin (TE19) with data mean (TE31)");
    hTr19 = fAnalysisMc->CreateTrigEff(var, 19);
    TArrayD tot31 = fAnalysisData->TotalTrigEff(31);
    for(Int_t bin = 1; bin <= hTr19->GetNbinsX(); ++bin){
      hTr19->SetBinContent(bin, hTr19->GetBinContent(bin) * tot31[0]);
      hTr19->SetBinError(bin, 0.);
    }
    hTr19->SetTitle("#epsilon_{trig}: 19 from MC, 31 mean from data");
    fAnalysisData->PrintTable(hTr19, TString("trigEff from MC (TE19) & data (TE31) ")+=var);
  }
  return hTr19;
}

//____________________________________
 TArrayD GFDstarXSections::TotalTrigEff() const
{
  TArrayD trigEff(2);
  if(this->GetTrigger() == 83){
    trigEff = this->EffTrigFromMc() ? 
      fAnalysisMc->TotalTrigEff() : fAnalysisData->TotalTrigEff();
  } else if(this->GetTrigger() == 84){
    this->Warning("TotalTrigEff", "Mixing MC (TE19) with data (TE31)");

    TArrayD mc19 = fAnalysisMc->TotalTrigEff(19); 
    TArrayD data31 = fAnalysisData->TotalTrigEff(31);
    trigEff[0] =  mc19[0] * data31[0];
    trigEff[1] = 0.;
    fAnalysisData->PrintPlusMinus(trigEff, "Total trigEff from MC (TE19) & data (TE31)");
  }
  return trigEff;
}


//____________________________________
 TH1* GFDstarXSections::DataVsMcDstar(const char* var, Bool_t perBinW, Bool_t log)
{
  // "pt", "eta, phi
  // returns data hist
  const Double_t allData = fAnalysisData->DrawDm()[0];
  TH1* histData = fAnalysisData->DrawDm(var, perBinW);
  const Double_t allMc = fAnalysisMc->DrawDm()[0];
  TH1* histMc = fAnalysisMc->DrawDm(var, perBinW);
  TH1* histMc2 = (fAnalysisOtherMc ? fAnalysisOtherMc->DrawDm(var, perBinW) : NULL);
  const Double_t allMc2 = (fAnalysisOtherMc ? fAnalysisOtherMc->DrawDm()[0] : 0.);

  if(!histMc || !histData){
    this->Error("DataVsMcDstar", "One hist missing!");
    return NULL;
  }

  GFHistManip::Scale(histData, 1./allData);
  GFHistManip::Scale(histMc, 1./allMc);
  if(histMc2) {
    GFHistManip::Scale(histMc2, 1./allMc2);
  }

  return this->DrawDataVsMC(histData, histMc, histMc2, perBinW, log);
}

//____________________________________
 TH1* GFDstarXSections::DataVsMcDstar(const char* var1, const char* var2, 
				     Int_t firstBin, Int_t lastBin, Bool_t perBinW, Bool_t log)
{
  // pt, eta, wGammaP
  const Double_t allData = fAnalysisData->DrawDm()[0];
  TH1* histData = fAnalysisData->DrawDm(var1, var2, firstBin, lastBin, perBinW);
  const Double_t allMc = fAnalysisMc->DrawDm()[0];
  TH1* histMc = fAnalysisMc->DrawDm(var1, var2, firstBin, lastBin, perBinW);
  TH1* histMc2 = (fAnalysisOtherMc ? 
		  fAnalysisOtherMc->DrawDm(var1, var2, firstBin, lastBin, perBinW) : NULL);
  const Double_t allMc2 = (fAnalysisOtherMc ? fAnalysisOtherMc->DrawDm()[0] : 0.);

  if(!histMc || !histData){
    this->Error("DataVsMcDstar", "One hist missing!");
    return NULL;
  }

  GFHistManip::Scale(histData, 1./allData);
  GFHistManip::Scale(histMc, 1./allMc);
  if(histMc2) {
    GFHistManip::Scale(histMc2, 1./allMc2);
  }

  return this->DrawDataVsMC(histData, histMc, histMc2, perBinW, log);
}

//____________________________________
 TH1* GFDstarXSections::DataVsMcDsDiJet(const char* var, Bool_t perBinW, Bool_t log)
{
  // Pt, Eta, Phi + Ds, Jet; xGam, Dphi, Deta,...
  const UInt_t dmRefFlag = 2; // inclusive data as reference
  const Double_t allData = fAnalysisData->DrawDmDsDiJet(dmRefFlag).At(0);
  TH1* histData = fAnalysisData->DrawDmDsDiJet(var, perBinW, dmRefFlag);
  const Double_t allMc   = fAnalysisMc->  DrawDmDsDiJet(dmRefFlag).At(0);
  TH1* histMc = fAnalysisMc->DrawDmDsDiJet(var, perBinW, dmRefFlag);
  const Double_t allMc2 = (fAnalysisOtherMc ? 
			   fAnalysisOtherMc->DrawDmDsDiJet(dmRefFlag).At(0) : 0.);
  TH1* histMc2 = (fAnalysisOtherMc ? 
		  fAnalysisOtherMc->DrawDmDsDiJet(var, perBinW, dmRefFlag) : NULL);

  if(!histMc || !histData){
    this->Error("DataVsMcDsDiJet", "One hist missing!");
    return NULL;
  }

  GFHistManip::Scale(histData, 1./allData);
  GFHistManip::Scale(histMc, 1./allMc);
  if(histMc2) {
    GFHistManip::Scale(histMc2, 1./allMc2);
  }

  TH1* hResult = this->DrawDataVsMC(histData, histMc, histMc2, perBinW, log);
  TString yTitle = hResult->GetYaxis()->GetTitle();
  Ssiz_t ind = yTitle.Index("D*");
  if(ind != kNPOS) yTitle.Insert(ind+2, " dijet");
  hResult->SetYTitle(yTitle);

  return hResult;
}

//____________________________________
 TH1* GFDstarXSections::DataVsMcDsJet(const char* var, const char* forwBack,
				     Bool_t perBinW, Bool_t log)
{
  // Pt, Eta, Phi + Ds, Jet; xGam, Dphi, Deta,...
  const UInt_t dmRefFlag = 2; // inclusive data as reference
  const Double_t allData = fAnalysisData->DrawDmDsJet(forwBack, dmRefFlag).At(0);
  TH1* histData = fAnalysisData->DrawDmDs1Jet(var, forwBack, perBinW, dmRefFlag);
  const Double_t allMc   = fAnalysisMc->  DrawDmDsJet(forwBack, dmRefFlag).At(0);
  TH1* histMc = fAnalysisMc->DrawDmDs1Jet(var, forwBack, perBinW, dmRefFlag);
  const Double_t allMc2 = (fAnalysisOtherMc ? 
			   fAnalysisOtherMc->DrawDmDsJet(forwBack, dmRefFlag).At(0) : 0.);
  TH1* histMc2 = (fAnalysisOtherMc ? 
		  fAnalysisOtherMc->DrawDmDs1Jet(var, forwBack, perBinW, dmRefFlag) : NULL);

  if(!histMc || !histData){
    this->Error("DataVsMcDsJet", "One hist missing!");
    return NULL;
  }

  GFHistManip::Scale(histData, 1./allData);
  GFHistManip::Scale(histMc, 1./allMc);
  if(histMc2) {
    GFHistManip::Scale(histMc2, 1./allMc2);
  }

  TH1* hResult = this->DrawDataVsMC(histData, histMc, histMc2, perBinW, log);
  TString yTitle = hResult->GetYaxis()->GetTitle();
  Ssiz_t ind = yTitle.Index("D*");
  if(ind != kNPOS) yTitle.Insert(ind+2, "+jet");
  hResult->SetYTitle(yTitle);

  return hResult;
}

//____________________________________
 TH1* GFDstarXSections::DrawDataVsMC(TH1* histData, TH1* histMc, TH1* histMc2,
				    Bool_t perBinW, Bool_t log)
{
  if(!histData || !histMc) return NULL;

  fHistManager->Clear();

  fHistManager->AddHist(histData,  1, fAnalysisData->GetPeriod(), "lp");
  histData->SetYTitle(perBinW ? "(n(D*)/N)/bin" : "n(D*)/N");
  histMc->SetOption("E2");
  fHistManager->AddHistSame(histMc, 1, 0, fAnalysisMc->GetPeriod(),"f");
  if(histMc2) {
    histMc2->SetOption("E2");
    fHistManager->AddHistSame(histMc2, 1, 0, fAnalysisOtherMc->GetPeriod(),"f");
  }
  //  fHistManager->AddLegend(0,0,fAnalysisData->GetSampleName());

  fHistManager->AddHistSame(histMc,  1, 0); // draw on top again
  fHistManager->AddHistSame(histData,  1, 0); // draw on top again

  if(log) fHistManager->SetLogY(1);
  //  TH1* hSb = this->GetDstarHistsAnalysisData()->GetHistManager()->GetHistsOf(2, 1)->At(0);
  if(this->GetDstarHistsAnalysisData()->GetHistManager()->GetHistsOf(2, 1)){
    fHistManager->AddHist(this->GetDstarHistsAnalysisData()->GetHistManager()
			  ->GetHistsOf(2, 1)->At(0), 0);
  }
  
  fHistManager->Draw();

  histData->SetLineStyle(1);
  histData->SetLineColor(kBlack);
  histData->SetMarkerColor(histData->GetLineColor());
  histData->SetMarkerStyle(kFullDotLarge);
  //  histData->SetMarkerSize(1.3);

  histMc->SetFillColor(8); // 8: dark green 5: yellow
  histMc->SetLineColor(histMc->GetFillColor());
  histMc->SetMarkerColor(histMc->GetFillColor());
  histMc->SetLineStyle(1);

  if(histMc2){
    histMc2->SetFillColor(kMagenta); // 8: dark green 5: yellow
    histMc2->SetLineColor(histMc2->GetFillColor());
    histMc2->SetMarkerColor(histMc2->GetFillColor());
    histMc2->SetLineStyle(1);
  }

  return histData;
}


//____________________________________
 void GFDstarXSections::DataVsMcTracks(const char* var)
{
  // Pt,Theta,Length,Phi,LH, etc.
  TString v(var);
  fHistMc1 = fAnalysisMc->DrawTrack(v+"K");
  fHistMc2 = fAnalysisMc->DrawTrack(v+"Pi");
  fHistMc3 = fAnalysisMc->DrawTrack(v+"Pis");

  fHistData1 = fAnalysisData->DrawTrack(v+"K");
  fHistData2 = fAnalysisData->DrawTrack(v+"Pi");
  fHistData3 = fAnalysisData->DrawTrack(v+"Pis");

  this->DataVsMcTracks();
}

//____________________________________
 void GFDstarXSections::DataVsMcTracks()
{
  GFHistManip::Normalise1stTo2nd(fHistMc1, fHistData1);//, "width"?
  GFHistManip::Normalise1stTo2nd(fHistMc2, fHistData2);
  GFHistManip::Normalise1stTo2nd(fHistMc3, fHistData3);

  fHistManager->Clear();
  fHistManager->AddHist(fHistMc1,0,fAnalysisMc->GetPeriod(), "f");
  fHistManager->AddLegend(0,0,fAnalysisMc->GetSampleName());
  fHistManager->AddHist(fHistMc2,0,fAnalysisMc->GetPeriod(), "f");
  fHistManager->AddLegend(0,1,fAnalysisMc->GetSampleName());
  fHistManager->AddHist(fHistMc3,0,fAnalysisMc->GetPeriod(), "f");
  fHistManager->AddLegend(0,2,fAnalysisMc->GetSampleName());
  fHistManager->SetHistsOption("E2");
  fHistManager->SetHistsFillColor(5); // 5: yellow

  fHistManager->AddHistSame(fHistData1,0,0,fAnalysisData->GetPeriod());
  fHistManager->AddHistSame(fHistData2,0,1,fAnalysisData->GetPeriod());
  fHistManager->AddHistSame(fHistData3,0,2,fAnalysisData->GetPeriod());

  fHistManager->SetNumHistsX(1);
  fHistManager->SetNumHistsY(3);
  fHistManager->Draw();

  this->SetHistPointerNull();
}

//____________________________________
 TH1* GFDstarXSections::DataVsMcTrack(const char* varTrack)
{
  // ptK, LHPi, thetaPis,...

  TString var(varTrack);
  const Double_t allData = fAnalysisData->DrawDm()[0];
  const Double_t allMc = fAnalysisMc->DrawDm()[0];
  const Double_t allMc2 = fAnalysisOtherMc ? fAnalysisOtherMc->DrawDm()[0] : 0.;

  TH1* histMc = fAnalysisMc->DrawTrack(var);
  TH1* histData = fAnalysisData->DrawTrack(var);
  TH1* histMc2 = (fAnalysisOtherMc ? fAnalysisOtherMc->DrawTrack(var) : NULL);

  if(!histMc || ! histData){
    this->Error("DataVsMcTrack", "One hist missing!");
    return NULL;
  }

  GFHistManip::Scale(histData, 1./allData);
  GFHistManip::Scale(histMc, 1./allMc);
  if(histMc2) {
    GFHistManip::Scale(histMc2, 1./allMc2);
  }

  return this->DrawDataVsMC(histData, histMc, histMc2, kFALSE, kFALSE);// perBinW,log
}

//____________________________________
 void GFDstarXSections::DataVsMcPtDsOverSumEt()
{
  fHistMc1 = fAnalysisMc->DrawDmPtDsOverSumEt();
  fHistData1 = fAnalysisData->DrawDmPtDsOverSumEt();

  GFHistManip::Normalise1stTo2nd(fHistMc1, fHistData1);//, "width"?

  fHistManager->Clear();
  fHistManager->AddHist(fHistMc1,0);
  fHistManager->SetHistsOption("E2");
  fHistManager->SetHistsFillColor(5); // 5: yellow
  fHistManager->AddHistSame(fHistData1,0,0);
  fHistManager->Draw();

  this->SetHistPointerNull();
}


//____________________________________
 TH1* GFDstarXSections::DataVsMcTrigEff(const char* var, Int_t TE)
{
  //pt,eta,phi,wGammaP; TE: 19/31/default=-1
  // retunrs the hist defining the axis limits

  TH1* histData = fAnalysisData->CreateTrigEff(var, TE);
  TGraphAsymmErrors* errorsData = fAnalysisData->TmpGraph();

  TH1* histAllMc = fAnalysisMc->CreateTrigEff(var, TE);
  TGraphAsymmErrors* errorsMc = fAnalysisMc->TmpGraph();

  if(!histData || !histAllMc) return NULL;

//   GFHistArray mc;
//   mc.AddAll(fAnalysisMc->GetHistsTmp());
//   fHistManager->AddHist(histData, 0,  fAnalysisData->GetPeriod());
//   for(Int_t i = 0; i < mc.GetEntriesFast(); ++i){
//     fHistManager->AddHistSame(mc[i], 0, 0, fAnalysisMc->GetSingleMc(i)->GetPeriod());
//   }
//   if(mc.GetEntriesFast() != (Int_t) fAnalysisMc->GetNumMc()){
//     this->Warning("DataVsMcTrigEff", "number of single hists not matching");
//   }
  TH1* hLowerMC = GFHistManip::CreateAnalogHist(histAllMc, errorsMc, kFALSE);
  TH1* hUpperMC = GFHistManip::CreateAnalogHist(histAllMc, errorsMc, kTRUE);
  TString name(Form("effTrigMc%s%sS%d", fAnalysisMc->GetPeriod(), var, this->GetTrigger()));
  GFHistManip::MakeUniqueName(name);
  TH1* histAllMcClone = GFHistManip::CreateAnalog(histAllMc, name, kFALSE, kTRUE);
  for(Int_t i = 0; i <= histAllMc->GetNbinsX()+1; ++i){ // under and overflow, too
    if(hLowerMC && hUpperMC) {
      histAllMc->SetBinError(i, FLT_MIN); // protect if no Graph
    }
    histData->SetBinError(i, FLT_MIN);// assuming same number of bins...
  }

  fHistManager->Clear();

  if(hLowerMC && hUpperMC) {
    fHistManager->AddHist(hUpperMC, 1);
    fHistManager->AddHistSame(histAllMc, 1, 0);
  } else {
    fHistManager->AddHist(histAllMc, 1);
    histAllMc->SetDrawOption("E");
  }
  fHistManager->AddHistSame(hLowerMC, 1, 0);
  fHistManager->AddHistSame(histData, 1,  0);
  fHistManager->AddObject(errorsData, 1, 0, "P");
  TString header;
  switch(TE){
  case 19:
    header = "DCRPh_Tc";
    break;
  case 31:
    header = "zVtx_sig_1";
    break;
  case -1:
    header = "DCRPh_Tc && zVtx_sig_1";
  }

  TLegend* leg = fHistManager->AddLegend(1, 0, header.Data(), kFALSE);
  leg->AddEntry(errorsData, fAnalysisData->GetPeriod(), "LP");
  leg->AddEntry(histAllMcClone, fAnalysisMc->GetPeriod(),"LF");

  if(hLowerMC && hUpperMC) {
    hUpperMC->SetFillStyle(1001);
    hUpperMC->SetFillColor(kYellow);
    hUpperMC->SetLineColor(10);
    hUpperMC->SetLineWidth(1);
 
    hLowerMC->SetFillStyle(1001);
    hLowerMC->SetFillColor(10);
    hLowerMC->SetLineColor(10);
    hLowerMC->SetLineWidth(1);
    hLowerMC->SetOption("HIST"); // some how needed   ????
  }

  histAllMc->SetMarkerColor(kRed);
  histAllMc->SetLineColor(kRed);
  histAllMc->SetLineStyle(2);

  histAllMcClone->SetFillColor(kYellow);
  histAllMc->TAttLine::Copy(*histAllMcClone);
  //  if(hLowerMC && hUpperMC) hUpperMC->TAttFill::Copy(*histAllMcClone);

  errorsData->SetMarkerStyle(8);
  errorsData->SetMarkerSize(1.3);

  const Bool_t oldDiff = fHistManager->DrawDiffStyle(kFALSE);
  fHistManager->Draw();
  fHistManager->DrawDiffStyle(oldDiff);

  if(hLowerMC && hUpperMC) {
    GFHistArray hists;
    hists.Add(hUpperMC); hists.Add(histData); hists.Add(hLowerMC);
    hUpperMC->SetMaximum(TMath::Max(1.01, GFHistManip::MinMaxOfHists(&hists, 1)));
    hUpperMC->SetMinimum(TMath::Min(.85, GFHistManip::MinMaxOfHists(&hists, -1)));
  }

  if(gPad) gPad->RedrawAxis();

  if(hLowerMC && hUpperMC) return hUpperMC;
  else return histAllMc;
}

//____________________________________
 void GFDstarXSections::DataVsMcL4EffSR(const char* var, Int_t dstar, Int_t fitModeData)
{
  // comparing L4 eff in signal region
  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){
    TH1* histData = fAnalysisData->CreateL4EffSR(vars[iVar]->GetName(), dstar, fitModeData);
    TGraphAsymmErrors* grData = fAnalysisData->TmpGraph();
    if(grData){
      histData->SetMarkerStyle(8);
      histData->TAttMarker::Copy(*grData);
      histData->TAttLine::Copy(*grData);
    }
    TH1* histAllMc = fAnalysisMc->CreateL4EffSR(vars[iVar]->GetName(), dstar);

    fHistManager->AddHist(histData, 0,  fAnalysisData->GetPeriod(), "LP");
    if(grData) fHistManager->AddObject(grData, 0, iVar, "P"); 
    fHistManager->AddHistSame(histAllMc, 0, iVar,  fAnalysisMc->GetPeriod());

    fHistManager->AddHist(static_cast<TH1*>(histAllMc->Clone()), 1, "All MC");
//     TLegend* leg = fHistManager->AddLegend(0, 0, NULL, kFALSE);
//     leg->AddEntry(grData, fAnalysisData->GetPeriod(), "LP"); 
    
//     for(Int_t j = 0; j < mc.GetEntriesFast(); ++j){ // before 'per lumi'
//       fHistManager->AddHistSame(mc[j], 1, iVar, fAnalysisMc->GetSingleMc(j)->GetPeriod());
//     }
  }
  
  fHistManager->SetHistsMinMax(0.8, 1.02);
  fHistManager->Draw();
  vars.Delete();
}


//____________________________________
 void GFDstarXSections::SetTrigger(Int_t trigger) 
{
  fTrigger = trigger;
  fAnalysisData->SetTrigger(fTrigger);
  fAnalysisMc->SetTrigger(fTrigger);
  if(fAnalysisOtherMc) fAnalysisOtherMc->SetTrigger(fTrigger);
  fInclXsec[0] = fInclXsec[1] = 0.; // invalidate
}

//____________________________________
 void GFDstarXSections::SetBatch(Bool_t yesNo)
{
  fHistManager->SetBatch(yesNo);

  fAnalysisData->SetBatch(yesNo);
  fAnalysisMc->SetBatch(yesNo);
  if(fAnalysisOtherMc) fAnalysisOtherMc->SetBatch(yesNo);
}

//____________________________________
 void GFDstarXSections::SetHistPointerNull()
{
  fHistNumDstar = fHistTrigEff = fHistRecEff = fHistAccept 
    = NULL;
  fHistData1 = fHistData2 = fHistData3 = fHistMc1 = fHistMc2 = fHistMc3
    = NULL;
}

//____________________________________
 Double_t GFDstarXSections::NormaliseToData(TH1* histMc, GFDstarHistsAnalysisMc* mc,
					   Bool_t acceptOnMc) const
{
  // returns scale that is used
  if(!histMc || !mc) {
    this->Error("NormaliseToData", "hist or GFDstarHistsAnalysisMc are missing");
    return 0.;
  }
  Double_t xSectionMc = mc->XSectionTot(acceptOnMc).At(0);
  Double_t scale = xSectionMc ? this->XSectionTot(acceptOnMc).At(0) / xSectionMc : 0.;
  histMc->Scale(scale);
  if(histMc->GetMaximumStored() != -1111.) histMc->SetMaximum(histMc->GetMaximum() * scale);
  return scale;
}

//____________________________________
 GFHistArray* GFDstarXSections::FitFractions(TH1* histXSection, GFHistArray* mcHists)
{
  // fitting the fraction of which fit best to sum up 'mcHists' to 'histXSection'
  // creates (!) output array with histXSection and clones (!) of mcHists!
  GFHistArray* result = new GFHistArray;
  result->Add(histXSection);

  TArrayD fractions = fHistFracFitter.Fit(histXSection, mcHists);

  TArrayD bins(*(histXSection->GetXaxis()->GetXbins())); // FIXME: if bins not variable!

  TH1* sumMc = new TH1F("sumMcFitted", "AllMc", bins.GetSize()-1, bins.GetArray());
  sumMc->SetOption("HIST");
  sumMc->SetLineWidth(2);
  result->Add(sumMc);

  for(Int_t i = 0; i < mcHists->GetEntriesFast(); ++i){
    TH1* h = static_cast<TH1*>((*mcHists)[i]->Clone());
    h->Scale(fractions[i]);
    result->Add(h);
    sumMc->Add(h);
  }

  return result;
}

//____________________________________
 void  GFDstarXSections::SetMergeNBins(Int_t n)
{
  fAnalysisData->GetMultiDmFitter()->SetMergeNBins(n);
  fAnalysisMc->GetMultiDmFitter()->SetMergeNBins(n);
  if(fAnalysisOtherMc) fAnalysisOtherMc->GetMultiDmFitter()->SetMergeNBins(n);
}

//____________________________________
 void  GFDstarXSections::SetLumi(const char* per)
{
  // set lumi array: 99ePlus, 2000 or All are known periods
  // if no period is given try to extract from file name of data
  fLumi.Set(2); // [1]: S83, [0]: S84

  TString period(per);
  if(!per){
    const TFile* dataFile = fAnalysisData->GetInFile();
    if(dataFile) period = dataFile->GetName();
  }

  if(period.Contains("99ePlus", TString::kIgnoreCase)) period = "99ePlus";
  if(period.Contains("2000")) period = "2000";
  if(period.Contains("All", TString::kIgnoreCase)) period = "";

  GFDstarRunPeriods p(83, period);
  fLumi[1] = p.GetTotalLumi();
  p.SetTrigger(84, period);
  fLumi[0] = p.GetTotalLumi();
}


//____________________________________
 void  GFDstarXSections::SetFitOrNotFlag(Int_t flag)
{
  fAnalysisData->GetMultiDmFitter()->SetFitOrNotFlag(flag);
  fAnalysisMc->GetMultiDmFitter()->SetFitOrNotFlag(flag);
  if(fAnalysisOtherMc) fAnalysisOtherMc->GetMultiDmFitter()->SetFitOrNotFlag(flag);
}

//____________________________________
 const char* GFDstarXSections::GetUnit(const TString &titleShort) const
{
  if(titleShort.Contains("[GeV]")) return "GeV";
  if(titleShort.Contains("[#circ]")) return "#circ";

  return NULL;
}