Plot/src/PlotTask.cc

// $Id: PlotTask.cc,v 1.1.2.2 2010/10/12 21:25:03 paus Exp $

#include <vector>
#include <TROOT.h>
#include <TSystem.h>
#include <TMath.h>
#include <TPad.h>
#include <TFile.h>
#include <TH1.h>
#include <TH1D.h>
#include <TBox.h>
#include <TMarker.h>
#include <TLatex.h>
#include "MitAna/DataUtil/interface/Debug.h"
#include "MitPlots/Style/interface/MitStyle.h"
#include "MitPlots/Plot/interface/PlotTask.h"

ClassImp(mithep::PlotTask)

using namespace std;
using namespace mithep;

//--------------------------------------------------------------------------------------------------
PlotTask::PlotTask(const TaskSamples *taskSamples, const double lumi) :
  fTask        (taskSamples),
  fHistStyles  (0),
  fTargetLumi  (lumi),
  fIdxHistMax  (-1),
  fNRebin      (1),
  fEmptyHist   (0),
  fDataHist    (0),
  fHistMinimum (0),
  fHistMaximum (0),
  fHistXMinimum(0),
  fHistXMaximum(0),
  fAxisTitleX  (""),
  fAxisTitleY  ("Number of Events"),
  fXLegend     (70.),
  fYLegend     (98.)
{
  // Constructor
}

//--------------------------------------------------------------------------------------------------
PlotTask::~PlotTask()
{
  // Destructor

  if (fEmptyHist)
    delete fEmptyHist;
  if (fDataHist)
    delete fDataHist;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::SetAxisTitles(const char* xtit, const char* ytit)
{
  fAxisTitleX = TString(xtit);
  fAxisTitleY = TString(ytit);

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::PlotContributions(const char* dir, const char* hist)
{
  // Show present list of defined samples

  // scale the histograms
  ScaleHistograms(dir,hist);
  FindHistMaximum();

  // say what we are doing
  printf("\n ==== Plotting Contributions -- %s ====\n\n",fTask->Name()->Data());
  printf("  index of highest histogram %d (0-%d)\n\n",fIdxHistMax,fHists.size()-1);

  // check x-axis title
  if (fAxisTitleX == TString(""))
    fAxisTitleX = TString(hist);

  // initialize the largest histogram properly and draw it
  TH1D* hmax = 0;
  if (fIdxHistMax == fHists.size())
    hmax = fDataHist;
  else
    hmax = fHists[fHists.size()-1];
  MitStyle::InitHist(hmax,fAxisTitleX.Data(),fAxisTitleY.Data(),kBlack);
  hmax->GetXaxis()->SetNdivisions(505);
  if (fHistMinimum != 0)
    hmax->SetMinimum(fHistMinimum);
  if (fIdxHistMax == fHists.size())
    hmax->Draw("E");
  else
    hmax->Draw("hist");

  // loop through samples and draw all histograms
  int iCol  = (EColor) kBlack;
  int iFill = 3001;
  for (UInt_t i=0; i<*fTask->NSamples(); i++) {
    //const Sample *s = fTask->GetSample(i);
    MitStyle::InitHist(fHists[i],hist,"Number of Events",(EColor) iCol);

    fHists[i]->SetLineColor(iCol);
    //fHists[i]->SetFillColor(iCol);
    //fHists[i]->SetFillStyle(iFill);
    fHists[i]->Draw("same;Hist");

    iFill++; iCol++;
  }

  // overlay a frame to put the text and boxes on
  OverlayFrame();

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::PlotStack(const char* dir, const char* hist)
{
  // Show present list of defined samples

  // scale the histograms
  ScaleHistograms(dir,hist);
  FindHistMaximum();

  // ensure the histogram styles are ready
  if (! fHistStyles)
    fHistStyles = new HistStyles();

  // check basics
  if (fStackedHists.size() < 1) {
    printf(" WARNING -- no histograms with name: %s. EXIT!\n",hist);
    return;
  }

  // determine the width of one bin
  double binW = fEmptyHist->GetBinWidth(1);
  char binWStr[7];
  sprintf(binWStr," /%7.3f",binW);

  // say what we are doing
  printf("\n ==== Plotting Stack -- %s ====\n",fTask->Name()->Data());
  printf("\n   bin width: %7.3f  in  xMin: %7.3f  <--> xMax: %7.3f\n\n",
         binW,fHistXMinimum,fHistXMaximum);

  // check x-axis title
  if (fAxisTitleX == TString(""))
    fAxisTitleX = TString(hist);

  // create the requested frame and initialize the empty histogram properly and draw it
  TH1D *hFrame = 0;
  if (fHistXMinimum != 0. || fHistXMaximum != 0.) {
    hFrame = new TH1D("Frame",fEmptyHist->GetTitle(),1,fHistXMinimum,fHistXMaximum);
  }
  else
    hFrame = fEmptyHist;

  MitStyle::InitHist(hFrame,fAxisTitleX.Data(),(fAxisTitleY+TString(binWStr)).Data(),kBlack);
  hFrame->GetXaxis()->SetNdivisions(505);
  if (fHistMinimum != 0)
    hFrame->SetMinimum(fHistMinimum);
  hFrame->SetMaximum(fHistMaximum*1.1);
  hFrame->Draw("");

  // loop through samples and draw all histograms
  fHistStyles->ResetStyle();
  for (UInt_t i=fStackedHists.size(); i>0; i--) {
    const Sample *s = fTask->GetSample(i);
    MitStyle::InitHist(fStackedHists[i-1],fAxisTitleX.Data(),fAxisTitleY.Data());
    if (i == fStackedHists.size()) {
      fHistStyles->ApplyCurrentStyle(fStackedHists[i-1]);
      fStackedHists[i-1]->Draw("same;Hist");
      fHistStyles->NextStyle();
    }
    else {
      if (*s->Legend() != TString(" ")) {
        // white out the histogram first
        fStackedHists[i-1]->SetFillColor(10);
        fStackedHists[i-1]->SetFillStyle(1001);
        fStackedHists[i-1]->DrawCopy("same;Hist");
        // overlay the histogram with proper hatching
        fHistStyles->ApplyCurrentStyle(fStackedHists[i-1]);
        fStackedHists[i-1]->Draw("same;Hist");
        fHistStyles->NextStyle();
      }
    }
  }

  if (fDataHist) {
    fHistStyles->ApplyDataStyle(fDataHist);
    fDataHist->Draw("same;E");
  }

  // overlay a frame to put the text and boxes on
  OverlayFrame();

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::ScaleHistograms(const char* dir, const char* hist)
{
  // Scale the histograms according to the cross section and the desired lumi and store them
  // for later use

  if (fHists.size() > 0) {
    printf(" WARNING - scaled histograms already exist. EXIT.");
    return;
  }

  fIdxHistMax = -1;

  // some useful definitions
  TString *dirFwk  = new TString("AnaFwkMod");
  TString *allEvts = new TString("hDAllEvents");
  TString  slash  ("/");

  // say what we are doing
  printf("\n ==== Extracting and Scaling Contributions -- %s ====\n\n",fTask->Name()->Data());
  printf("  target luminosity %f 1/fb\n\n",fTargetLumi/1000.);

  //------------------------------------------------------------------------------------------------
  // go through the Monte Carlo samples
  //------------------------------------------------------------------------------------------------
  // loop through samples and determine maximum
  printf("\n monte carlo \n");
  for (UInt_t i=0; i<*fTask->NSamples(); i++) {
    const Sample *s = fTask->GetSample(i);
    // open file belonging to this sample
    TFile *fif = new TFile((*fTask->Dir()+slash+*s->File()).Data());
    // make sure the file exists
    if (fif->IsOpen() == kFALSE) {
      printf(" WARNING -- sample  %s  does not have a histogram file. Continue without!\n",
             s->Name()->Data());
      continue;
    }
    // read and determine general properties of this sample
    TDirectory *dirTmp = (TDirectory*) gROOT->FindObject(dirFwk->Data());
    if (dirTmp)
      fif->cd(dirFwk->Data());
    TH1D *hAllEvts = (TH1D*) gROOT->FindObject(allEvts->Data());
    if (! hAllEvts) {
      printf(" WARNING -- sample  %s  does not have a framework file. Next sample!\n",
             s->Name()->Data());
      continue;
    }
    double nEvts   = hAllEvts->GetEntries();
    double lumi    = nEvts / *s->Xsec();
    double factor,scale = *s->Scale();
    if      (lumi < 0)
      factor = 1.0;
    else if (lumi > 0)
      factor = fTargetLumi/lumi;
    else
      factor = 0;

    TH1D *h    = (TH1D*) fif->FindObjectAny(hist);
    if (! h) {
      printf(" WARNING -- sample  %s  does not have requested histogram. Next sample!\n",
             s->Name()->Data());
      continue;
    }
    if (! fEmptyHist) {
      fEmptyHist = new TH1D(factor * scale * (*h));
      fEmptyHist->Rebin(fNRebin);
      fEmptyHist->Reset();
    }

    double nEvtsSelRaw = h->GetEntries();
    double nEvtsSel    = nEvtsSelRaw * factor * scale;
    double nEvtsSelErr = TMath::Sqrt(nEvtsSelRaw) * factor * scale;

    printf(" -> %-40s - %14.0f %10.3f +- %6.3f %16.7f: %16.4f (x %f x %f - %p)\n",
           s->Name()->Data(),nEvts,nEvtsSel,nEvtsSelErr,*s->Xsec(),lumi,factor,scale,(void*)h);

    // scale it
    TH1D *hTmp = new TH1D(factor * scale * (*h));
    hTmp->Rebin(fNRebin);
    // put it into our collection
    fHists.push_back(hTmp);
    // make a new histogram from the last added one and add the recent
    if (i == 0)
      fStackedHists.push_back(hTmp);
    else {
      TH1D *hStackedTmp = new TH1D(*fStackedHists[fStackedHists.size()-1]);
      // remember the stored histogram is already rebinned
      hStackedTmp->Add(hTmp);
      fStackedHists.push_back(hStackedTmp);
    }
  }

  //------------------------------------------------------------------------------------------------
  // go through the data samples
  //------------------------------------------------------------------------------------------------
  // loop through data samples and add them up, straight away
  if (*fTask->NDataSamples() > 0)
    printf("\n data \n");
  for (UInt_t i=0; i<*fTask->NDataSamples(); i++) {
    const Sample *s = fTask->GetDataSample(i);
    // open file belonging to the data sample
    TFile *fif = new TFile((*fTask->Dir()+slash+*s->File()).Data());
    // make sure the file exists
    if (fif->IsOpen() == kFALSE) {
      printf(" WARNING -- sample  %s  does not have a histogram file. Continue without!\n",
             s->Name()->Data());
      continue;
    }
    else {
      // read and determine general properties of this sample
      TDirectory *dirTmp = (TDirectory*) gROOT->FindObject(dirFwk->Data());
      if (dirTmp)
        fif->cd(dirFwk->Data());
      TH1D *hAllEvts = (TH1D*) gROOT->FindObject(allEvts->Data());
      if (! hAllEvts) {
        printf(" WARNING -- sample  %s  does not have a framework file. Next sample!\n",
               s->Name()->Data());
        continue;
      }
      else {
        double nEvts = hAllEvts->GetEntries();

        TH1D *h = (TH1D*) fif->FindObjectAny(hist);
        if (! h) {
          printf(" WARNING -- sample  %s  does not have requested histogram. Next sample!\n",
                 s->Name()->Data());
          continue;
        }
        else
          // rebin it
          h->Rebin(fNRebin);

        double nEvtsSelRaw = h->GetEntries();
        double nEvtsSel    = nEvtsSelRaw;
        double nEvtsSelErr = TMath::Sqrt(nEvtsSelRaw);

        printf(" -> %-40s - %14.0f %10.3f +- %6.3f %16.7f: %16.4f (x %f)\n",
               s->Name()->Data(),nEvts,nEvtsSel,nEvtsSelErr,*s->Xsec(),fTargetLumi,1.0);


        // construct the complete data histogram
        if (! fDataHist)
          fDataHist = new TH1D(1.0 * (*h));
        else
          fDataHist->Add(h);
      }
    }
  }

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::FindHistMaximum()
{
  // Find maximum of all histograms

  // first check whether value was overwritten by hand
  if (fHistMaximum>0.)
    return;
  
  // search through the single histograms
  for (UInt_t i=0; i<fHists.size(); i++) {
    if (fHists[i]->GetMaximum() > fHistMaximum) {
      fHistMaximum = fHists[i]->GetMaximum();
      fIdxHistMax  = i;
    }
  }

  // search through the stacked histograms
  fHistMaximum = 0.;
  for (UInt_t i=0; i<fStackedHists.size(); i++) {
    if (fStackedHists[i]->GetMaximum() > fHistMaximum)
      fHistMaximum = fStackedHists[i]->GetMaximum();
  }
  
  // check the data histogram if present
  if (fDataHist && fDataHist->GetMaximum() > fHistMaximum) {
    fHistMaximum = fDataHist->GetMaximum();
    fIdxHistMax  = fHists.size();
  }

  printf(" Histogram maximum is set to be: %f\n",fHistMaximum);

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::OverlayEmptyHist() const
{
  // Overlay an empty histogram onto the picture
  if (fEmptyHist)
    fEmptyHist->Draw("same");

  return;
}

//--------------------------------------------------------------------------------------------------
void PlotTask::OverlayFrame() const
{
  // Overlay a linear frame from user coordinates (0-100,0-100)

  // create new transparent pad for the text
  TPad *transPad = new TPad("transPad","Transparent Pad",0,0,1,1);
  transPad->SetFillStyle(4000);
  transPad->Draw();
  transPad->cd();
  // find out the right normalization to define the new range (histogram: 0,0 -> 100,100)
  double xtot = 1./(1. - transPad->GetLeftMargin() - transPad->GetRightMargin());
  double ytot = 1./(1. - transPad->GetBottomMargin() - transPad->GetTopMargin());
  transPad->Range((xtot*transPad->GetLeftMargin()  *-100      ),
                  (ytot*transPad->GetBottomMargin()*-100      ),
                  (xtot*transPad->GetRightMargin() * 100 + 100),
                  (ytot*transPad->GetTopMargin()   * 100 + 100));
  MDB(kGeneral,1)
    printf(" Range: %f %f %f %f\n",
           (xtot*transPad->GetLeftMargin()  *-100      ),
           (ytot*transPad->GetBottomMargin()*-100      ),
           (xtot*transPad->GetRightMargin() * 100 + 100),
           (ytot*transPad->GetTopMargin()   * 100 + 100) );

  // define the basic text to be used
  TLatex* text = new TLatex();
  text->SetTextFont (42);
  text->SetTextAlign(12);
  text->SetTextSize (0.03);

  // define the basic box to be used
  TBox *box = new TBox();
  box->SetLineWidth(2);
  // base coordinates for the text/boxes
  double xCorner = fXLegend, yCorner = fYLegend, yDelLine = 4.;
  double xIndent = 1.5*yDelLine;
  // count samples with non empty legend
  int nLegends = 0;
  for (UInt_t i=*fTask->NSamples(); i>0; i--) {
    // attach to the specific sample
    const Sample *s = fTask->GetSample(i-1);
    if (*s->Legend() != TString(" "))
      nLegends++;
  }
  int iDat = 0;
  if (fDataHist) {
    iDat = 1;
    nLegends++;
  }

  // set start values
  fHistStyles->ResetStyle();
  int iLeg  = 0;
  // loop through the sampels
  for (UInt_t i=*fTask->NSamples(); i>0; i--) {
    // attach to the specific sample
    const Sample *s = fTask->GetSample(i-1);
    // calculate corners for the text
    double xText = xCorner+xIndent, yText = yCorner-float((iLeg+iDat)+0.5)*yDelLine;
    // say what goes where
    MDB(kGeneral,1)
      printf(" Adding box at: (x1,y1,x2,y2) = (%6.2f,%6.2f,%6.2f,%6.2f)\n",
             xCorner+0.15*xIndent,yText-0.3*yDelLine,
             xCorner+0.85*xIndent,yText+0.4*yDelLine);
    // plot the box
    if (*s->Legend() != TString(" ")) {
      const HistStyle *hStyle = fHistStyles->CurrentStyle();
      box->SetFillStyle(0);
      box->SetFillColor(hStyle->Color());
      box->SetLineColor(hStyle->Color());
      box->DrawBox(xCorner+0.15*xIndent,yText-0.3*yDelLine,
                   xCorner+0.85*xIndent,yText+0.4*yDelLine);
      box->SetFillStyle(hStyle->FillStyle());
      box->DrawBox(xCorner+0.15*xIndent,yText-0.3*yDelLine,
                   xCorner+0.85*xIndent,yText+0.4*yDelLine);

      TString pText = *s->Legend();
      MDB(kGeneral,1)
        printf(" Adding text \"%s\" at: (x,y) = (%6.2f,%6.2f)\n",pText.Data(),xText,yText);
      // set the proper text color
      text->SetTextColor(hStyle->Color());
      // plot the text
      text->SetTextAlign(12);
      text->DrawLatex(xText,yText,pText.Data());

      fHistStyles->NextStyle();
      // keep track of how many were printed
      iLeg++;
    }
  }

  // attach to the data sample
  if (fDataHist) {
    const Sample *s = fTask->GetDataSample(0);
    // calculate corners for the text
    double xText = xCorner+xIndent, yText = yCorner-float(1-0.7)*yDelLine;
    // say what goes where
    MDB(kGeneral,1)
      printf(" Adding marker at: (x,y) = (%6.2f,%6.2f)\n",
             xCorner+0.5*xIndent,yText);
    // plot the marker
    if (*s->Legend() != TString(" ")) {
      const HistStyle *hStyle = fHistStyles->GetDataStyle();

      // draw marker
      TMarker *m = new TMarker(xCorner+0.5*xIndent,yText,20);
      m->SetMarkerColor(hStyle->Color      ());
      m->SetMarkerSize (hStyle->MarkerSize ());
      m->SetMarkerStyle(hStyle->MarkerStyle());
      m->Draw();
      // draw text
      char l[1024];
      sprintf(l,"%4.2f",fTargetLumi);
      //TString pText = *s->Legend() + TString(" (#int L = ") + TString(l) + TString("/pb)");
      TString pText = *s->Legend() + TString(" (") + TString(l) + TString("/pb)");
      MDB(kGeneral,1)
        printf(" Adding text \"%s\" at: (x,y) = (%6.2f,%6.2f)\n",pText.Data(),xText,yText);
      // set the proper text color
      text->SetTextColor(hStyle->Color());
      // plot the text
      text->SetTextAlign(12);
      text->DrawLatex(xText,yText,pText.Data());
    }
  }

  // Make sure the histogram frame is nice
  box->SetFillStyle(0);
  box->SetLineColor(kBlack);
  box->DrawBox(0,0,100,100);

  delete text;
  delete box;

  return;
}
Revision:	1.2
Committed:	Tue Jan 25 11:30:30 2011 UTC (14 years, 3 months ago) by paus
Content type:	text/plain
Branch:	MAIN
Changes since 1.1:	+532 -0 lines
Log Message:	2nd try
#	User	Rev	Content
1	paus	1.2	// $Id: PlotTask.cc,v 1.1.2.2 2010/10/12 21:25:03 paus Exp $
2
3			#include <vector>
4			#include <TROOT.h>
5			#include <TSystem.h>
6			#include <TMath.h>
7			#include <TPad.h>
8			#include <TFile.h>
9			#include <TH1.h>
10			#include <TH1D.h>
11			#include <TBox.h>
12			#include <TMarker.h>
13			#include <TLatex.h>
14			#include "MitAna/DataUtil/interface/Debug.h"
15			#include "MitPlots/Style/interface/MitStyle.h"
16			#include "MitPlots/Plot/interface/PlotTask.h"
17
18			ClassImp(mithep::PlotTask)
19
20			using namespace std;
21			using namespace mithep;
22
23			//--------------------------------------------------------------------------------------------------
24			PlotTask::PlotTask(const TaskSamples *taskSamples, const double lumi) :
25			fTask (taskSamples),
26			fHistStyles (0),
27			fTargetLumi (lumi),
28			fIdxHistMax (-1),
29			fNRebin (1),
30			fEmptyHist (0),
31			fDataHist (0),
32			fHistMinimum (0),
33			fHistMaximum (0),
34			fHistXMinimum(0),
35			fHistXMaximum(0),
36			fAxisTitleX (""),
37			fAxisTitleY ("Number of Events"),
38			fXLegend (70.),
39			fYLegend (98.)
40			{
41			// Constructor
42			}
43
44			//--------------------------------------------------------------------------------------------------
45			PlotTask::~PlotTask()
46			{
47			// Destructor
48
49			if (fEmptyHist)
50			delete fEmptyHist;
51			if (fDataHist)
52			delete fDataHist;
53			}
54
55			//--------------------------------------------------------------------------------------------------
56			void PlotTask::SetAxisTitles(const char* xtit, const char* ytit)
57			{
58			fAxisTitleX = TString(xtit);
59			fAxisTitleY = TString(ytit);
60
61			return;
62			}
63
64			//--------------------------------------------------------------------------------------------------
65			void PlotTask::PlotContributions(const char* dir, const char* hist)
66			{
67			// Show present list of defined samples
68
69			// scale the histograms
70			ScaleHistograms(dir,hist);
71			FindHistMaximum();
72
73			// say what we are doing
74			printf("\n ==== Plotting Contributions -- %s ====\n\n",fTask->Name()->Data());
75			printf(" index of highest histogram %d (0-%d)\n\n",fIdxHistMax,fHists.size()-1);
76
77			// check x-axis title
78			if (fAxisTitleX == TString(""))
79			fAxisTitleX = TString(hist);
80
81			// initialize the largest histogram properly and draw it
82			TH1D* hmax = 0;
83			if (fIdxHistMax == fHists.size())
84			hmax = fDataHist;
85			else
86			hmax = fHists[fHists.size()-1];
87			MitStyle::InitHist(hmax,fAxisTitleX.Data(),fAxisTitleY.Data(),kBlack);
88			hmax->GetXaxis()->SetNdivisions(505);
89			if (fHistMinimum != 0)
90			hmax->SetMinimum(fHistMinimum);
91			if (fIdxHistMax == fHists.size())
92			hmax->Draw("E");
93			else
94			hmax->Draw("hist");
95
96			// loop through samples and draw all histograms
97			int iCol = (EColor) kBlack;
98			int iFill = 3001;
99			for (UInt_t i=0; i<*fTask->NSamples(); i++) {
100			//const Sample *s = fTask->GetSample(i);
101			MitStyle::InitHist(fHists[i],hist,"Number of Events",(EColor) iCol);
102
103			fHists[i]->SetLineColor(iCol);
104			//fHists[i]->SetFillColor(iCol);
105			//fHists[i]->SetFillStyle(iFill);
106			fHists[i]->Draw("same;Hist");
107
108			iFill++; iCol++;
109			}
110
111			// overlay a frame to put the text and boxes on
112			OverlayFrame();
113
114			return;
115			}
116
117			//--------------------------------------------------------------------------------------------------
118			void PlotTask::PlotStack(const char* dir, const char* hist)
119			{
120			// Show present list of defined samples
121
122			// scale the histograms
123			ScaleHistograms(dir,hist);
124			FindHistMaximum();
125
126			// ensure the histogram styles are ready
127			if (! fHistStyles)
128			fHistStyles = new HistStyles();
129
130			// check basics
131			if (fStackedHists.size() < 1) {
132			printf(" WARNING -- no histograms with name: %s. EXIT!\n",hist);
133			return;
134			}
135
136			// determine the width of one bin
137			double binW = fEmptyHist->GetBinWidth(1);
138			char binWStr[7];
139			sprintf(binWStr," /%7.3f",binW);
140
141			// say what we are doing
142			printf("\n ==== Plotting Stack -- %s ====\n",fTask->Name()->Data());
143			printf("\n bin width: %7.3f in xMin: %7.3f <--> xMax: %7.3f\n\n",
144			binW,fHistXMinimum,fHistXMaximum);
145
146			// check x-axis title
147			if (fAxisTitleX == TString(""))
148			fAxisTitleX = TString(hist);
149
150			// create the requested frame and initialize the empty histogram properly and draw it
151			TH1D *hFrame = 0;
152			if (fHistXMinimum != 0. \|\| fHistXMaximum != 0.) {
153			hFrame = new TH1D("Frame",fEmptyHist->GetTitle(),1,fHistXMinimum,fHistXMaximum);
154			}
155			else
156			hFrame = fEmptyHist;
157
158			MitStyle::InitHist(hFrame,fAxisTitleX.Data(),(fAxisTitleY+TString(binWStr)).Data(),kBlack);
159			hFrame->GetXaxis()->SetNdivisions(505);
160			if (fHistMinimum != 0)
161			hFrame->SetMinimum(fHistMinimum);
162			hFrame->SetMaximum(fHistMaximum*1.1);
163			hFrame->Draw("");
164
165			// loop through samples and draw all histograms
166			fHistStyles->ResetStyle();
167			for (UInt_t i=fStackedHists.size(); i>0; i--) {
168			const Sample *s = fTask->GetSample(i);
169			MitStyle::InitHist(fStackedHists[i-1],fAxisTitleX.Data(),fAxisTitleY.Data());
170			if (i == fStackedHists.size()) {
171			fHistStyles->ApplyCurrentStyle(fStackedHists[i-1]);
172			fStackedHists[i-1]->Draw("same;Hist");
173			fHistStyles->NextStyle();
174			}
175			else {
176			if (*s->Legend() != TString(" ")) {
177			// white out the histogram first
178			fStackedHists[i-1]->SetFillColor(10);
179			fStackedHists[i-1]->SetFillStyle(1001);
180			fStackedHists[i-1]->DrawCopy("same;Hist");
181			// overlay the histogram with proper hatching
182			fHistStyles->ApplyCurrentStyle(fStackedHists[i-1]);
183			fStackedHists[i-1]->Draw("same;Hist");
184			fHistStyles->NextStyle();
185			}
186			}
187			}
188
189			if (fDataHist) {
190			fHistStyles->ApplyDataStyle(fDataHist);
191			fDataHist->Draw("same;E");
192			}
193
194			// overlay a frame to put the text and boxes on
195			OverlayFrame();
196
197			return;
198			}
199
200			//--------------------------------------------------------------------------------------------------
201			void PlotTask::ScaleHistograms(const char* dir, const char* hist)
202			{
203			// Scale the histograms according to the cross section and the desired lumi and store them
204			// for later use
205
206			if (fHists.size() > 0) {
207			printf(" WARNING - scaled histograms already exist. EXIT.");
208			return;
209			}
210
211			fIdxHistMax = -1;
212
213			// some useful definitions
214			TString *dirFwk = new TString("AnaFwkMod");
215			TString *allEvts = new TString("hDAllEvents");
216			TString slash ("/");
217
218			// say what we are doing
219			printf("\n ==== Extracting and Scaling Contributions -- %s ====\n\n",fTask->Name()->Data());
220			printf(" target luminosity %f 1/fb\n\n",fTargetLumi/1000.);
221
222			//------------------------------------------------------------------------------------------------
223			// go through the Monte Carlo samples
224			//------------------------------------------------------------------------------------------------
225			// loop through samples and determine maximum
226			printf("\n monte carlo \n");
227			for (UInt_t i=0; i<*fTask->NSamples(); i++) {
228			const Sample *s = fTask->GetSample(i);
229			// open file belonging to this sample
230			TFile fif = new TFile((fTask->Dir()+slash+*s->File()).Data());
231			// make sure the file exists
232			if (fif->IsOpen() == kFALSE) {
233			printf(" WARNING -- sample %s does not have a histogram file. Continue without!\n",
234			s->Name()->Data());
235			continue;
236			}
237			// read and determine general properties of this sample
238			TDirectory dirTmp = (TDirectory) gROOT->FindObject(dirFwk->Data());
239			if (dirTmp)
240			fif->cd(dirFwk->Data());
241			TH1D hAllEvts = (TH1D) gROOT->FindObject(allEvts->Data());
242			if (! hAllEvts) {
243			printf(" WARNING -- sample %s does not have a framework file. Next sample!\n",
244			s->Name()->Data());
245			continue;
246			}
247			double nEvts = hAllEvts->GetEntries();
248			double lumi = nEvts / *s->Xsec();
249			double factor,scale = *s->Scale();
250			if (lumi < 0)
251			factor = 1.0;
252			else if (lumi > 0)
253			factor = fTargetLumi/lumi;
254			else
255			factor = 0;
256
257			TH1D h = (TH1D) fif->FindObjectAny(hist);
258			if (! h) {
259			printf(" WARNING -- sample %s does not have requested histogram. Next sample!\n",
260			s->Name()->Data());
261			continue;
262			}
263			if (! fEmptyHist) {
264			fEmptyHist = new TH1D(factor * scale * (*h));
265			fEmptyHist->Rebin(fNRebin);
266			fEmptyHist->Reset();
267			}
268
269			double nEvtsSelRaw = h->GetEntries();
270			double nEvtsSel = nEvtsSelRaw * factor * scale;
271			double nEvtsSelErr = TMath::Sqrt(nEvtsSelRaw) * factor * scale;
272
273			printf(" -> %-40s - %14.0f %10.3f +- %6.3f %16.7f: %16.4f (x %f x %f - %p)\n",
274			s->Name()->Data(),nEvts,nEvtsSel,nEvtsSelErr,s->Xsec(),lumi,factor,scale,(void)h);
275
276			// scale it
277			TH1D hTmp = new TH1D(factor scale * (*h));
278			hTmp->Rebin(fNRebin);
279			// put it into our collection
280			fHists.push_back(hTmp);
281			// make a new histogram from the last added one and add the recent
282			if (i == 0)
283			fStackedHists.push_back(hTmp);
284			else {
285			TH1D hStackedTmp = new TH1D(fStackedHists[fStackedHists.size()-1]);
286			// remember the stored histogram is already rebinned
287			hStackedTmp->Add(hTmp);
288			fStackedHists.push_back(hStackedTmp);
289			}
290			}
291
292			//------------------------------------------------------------------------------------------------
293			// go through the data samples
294			//------------------------------------------------------------------------------------------------
295			// loop through data samples and add them up, straight away
296			if (*fTask->NDataSamples() > 0)
297			printf("\n data \n");
298			for (UInt_t i=0; i<*fTask->NDataSamples(); i++) {
299			const Sample *s = fTask->GetDataSample(i);
300			// open file belonging to the data sample
301			TFile fif = new TFile((fTask->Dir()+slash+*s->File()).Data());
302			// make sure the file exists
303			if (fif->IsOpen() == kFALSE) {
304			printf(" WARNING -- sample %s does not have a histogram file. Continue without!\n",
305			s->Name()->Data());
306			continue;
307			}
308			else {
309			// read and determine general properties of this sample
310			TDirectory dirTmp = (TDirectory) gROOT->FindObject(dirFwk->Data());
311			if (dirTmp)
312			fif->cd(dirFwk->Data());
313			TH1D hAllEvts = (TH1D) gROOT->FindObject(allEvts->Data());
314			if (! hAllEvts) {
315			printf(" WARNING -- sample %s does not have a framework file. Next sample!\n",
316			s->Name()->Data());
317			continue;
318			}
319			else {
320			double nEvts = hAllEvts->GetEntries();
321
322			TH1D h = (TH1D) fif->FindObjectAny(hist);
323			if (! h) {
324			printf(" WARNING -- sample %s does not have requested histogram. Next sample!\n",
325			s->Name()->Data());
326			continue;
327			}
328			else
329			// rebin it
330			h->Rebin(fNRebin);
331
332			double nEvtsSelRaw = h->GetEntries();
333			double nEvtsSel = nEvtsSelRaw;
334			double nEvtsSelErr = TMath::Sqrt(nEvtsSelRaw);
335
336			printf(" -> %-40s - %14.0f %10.3f +- %6.3f %16.7f: %16.4f (x %f)\n",
337			s->Name()->Data(),nEvts,nEvtsSel,nEvtsSelErr,*s->Xsec(),fTargetLumi,1.0);
338
339
340
341			// construct the complete data histogram
342			if (! fDataHist)
343			fDataHist = new TH1D(1.0 * (*h));
344			else
345			fDataHist->Add(h);
346			}
347			}
348			}
349
350			return;
351			}
352
353			//--------------------------------------------------------------------------------------------------
354			void PlotTask::FindHistMaximum()
355			{
356			// Find maximum of all histograms
357
358			// first check whether value was overwritten by hand
359			if (fHistMaximum>0.)
360			return;
361
362			// search through the single histograms
363			for (UInt_t i=0; i<fHists.size(); i++) {
364			if (fHists[i]->GetMaximum() > fHistMaximum) {
365			fHistMaximum = fHists[i]->GetMaximum();
366			fIdxHistMax = i;
367			}
368			}
369
370			// search through the stacked histograms
371			fHistMaximum = 0.;
372			for (UInt_t i=0; i<fStackedHists.size(); i++) {
373			if (fStackedHists[i]->GetMaximum() > fHistMaximum)
374			fHistMaximum = fStackedHists[i]->GetMaximum();
375			}
376
377			// check the data histogram if present
378			if (fDataHist && fDataHist->GetMaximum() > fHistMaximum) {
379			fHistMaximum = fDataHist->GetMaximum();
380			fIdxHistMax = fHists.size();
381			}
382
383			printf(" Histogram maximum is set to be: %f\n",fHistMaximum);
384
385			return;
386			}
387
388			//--------------------------------------------------------------------------------------------------
389			void PlotTask::OverlayEmptyHist() const
390			{
391			// Overlay an empty histogram onto the picture
392			if (fEmptyHist)
393			fEmptyHist->Draw("same");
394
395			return;
396			}
397
398			//--------------------------------------------------------------------------------------------------
399			void PlotTask::OverlayFrame() const
400			{
401			// Overlay a linear frame from user coordinates (0-100,0-100)
402
403			// create new transparent pad for the text
404			TPad *transPad = new TPad("transPad","Transparent Pad",0,0,1,1);
405			transPad->SetFillStyle(4000);
406			transPad->Draw();
407			transPad->cd();
408			// find out the right normalization to define the new range (histogram: 0,0 -> 100,100)
409			double xtot = 1./(1. - transPad->GetLeftMargin() - transPad->GetRightMargin());
410			double ytot = 1./(1. - transPad->GetBottomMargin() - transPad->GetTopMargin());
411			transPad->Range((xtottransPad->GetLeftMargin() -100 ),
412			(ytottransPad->GetBottomMargin()-100 ),
413			(xtottransPad->GetRightMargin() 100 + 100),
414			(ytottransPad->GetTopMargin() 100 + 100));
415			MDB(kGeneral,1)
416			printf(" Range: %f %f %f %f\n",
417			(xtottransPad->GetLeftMargin() -100 ),
418			(ytottransPad->GetBottomMargin()-100 ),
419			(xtottransPad->GetRightMargin() 100 + 100),
420			(ytottransPad->GetTopMargin() 100 + 100) );
421
422			// define the basic text to be used
423			TLatex* text = new TLatex();
424			text->SetTextFont (42);
425			text->SetTextAlign(12);
426			text->SetTextSize (0.03);
427
428			// define the basic box to be used
429			TBox *box = new TBox();
430			box->SetLineWidth(2);
431			// base coordinates for the text/boxes
432			double xCorner = fXLegend, yCorner = fYLegend, yDelLine = 4.;
433			double xIndent = 1.5*yDelLine;
434			// count samples with non empty legend
435			int nLegends = 0;
436			for (UInt_t i=*fTask->NSamples(); i>0; i--) {
437			// attach to the specific sample
438			const Sample *s = fTask->GetSample(i-1);
439			if (*s->Legend() != TString(" "))
440			nLegends++;
441			}
442			int iDat = 0;
443			if (fDataHist) {
444			iDat = 1;
445			nLegends++;
446			}
447
448			// set start values
449			fHistStyles->ResetStyle();
450			int iLeg = 0;
451			// loop through the sampels
452			for (UInt_t i=*fTask->NSamples(); i>0; i--) {
453			// attach to the specific sample
454			const Sample *s = fTask->GetSample(i-1);
455			// calculate corners for the text
456			double xText = xCorner+xIndent, yText = yCorner-float((iLeg+iDat)+0.5)*yDelLine;
457			// say what goes where
458			MDB(kGeneral,1)
459			printf(" Adding box at: (x1,y1,x2,y2) = (%6.2f,%6.2f,%6.2f,%6.2f)\n",
460			xCorner+0.15xIndent,yText-0.3yDelLine,
461			xCorner+0.85xIndent,yText+0.4yDelLine);
462			// plot the box
463			if (*s->Legend() != TString(" ")) {
464			const HistStyle *hStyle = fHistStyles->CurrentStyle();
465			box->SetFillStyle(0);
466			box->SetFillColor(hStyle->Color());
467			box->SetLineColor(hStyle->Color());
468			box->DrawBox(xCorner+0.15xIndent,yText-0.3yDelLine,
469			xCorner+0.85xIndent,yText+0.4yDelLine);
470			box->SetFillStyle(hStyle->FillStyle());
471			box->DrawBox(xCorner+0.15xIndent,yText-0.3yDelLine,
472			xCorner+0.85xIndent,yText+0.4yDelLine);
473
474			TString pText = *s->Legend();
475			MDB(kGeneral,1)
476			printf(" Adding text \"%s\" at: (x,y) = (%6.2f,%6.2f)\n",pText.Data(),xText,yText);
477			// set the proper text color
478			text->SetTextColor(hStyle->Color());
479			// plot the text
480			text->SetTextAlign(12);
481			text->DrawLatex(xText,yText,pText.Data());
482
483			fHistStyles->NextStyle();
484			// keep track of how many were printed
485			iLeg++;
486			}
487			}
488
489			// attach to the data sample
490			if (fDataHist) {
491			const Sample *s = fTask->GetDataSample(0);
492			// calculate corners for the text
493			double xText = xCorner+xIndent, yText = yCorner-float(1-0.7)*yDelLine;
494			// say what goes where
495			MDB(kGeneral,1)
496			printf(" Adding marker at: (x,y) = (%6.2f,%6.2f)\n",
497			xCorner+0.5*xIndent,yText);
498			// plot the marker
499			if (*s->Legend() != TString(" ")) {
500			const HistStyle *hStyle = fHistStyles->GetDataStyle();
501
502			// draw marker
503			TMarker m = new TMarker(xCorner+0.5xIndent,yText,20);
504			m->SetMarkerColor(hStyle->Color ());
505			m->SetMarkerSize (hStyle->MarkerSize ());
506			m->SetMarkerStyle(hStyle->MarkerStyle());
507			m->Draw();
508			// draw text
509			char l[1024];
510			sprintf(l,"%4.2f",fTargetLumi);
511			//TString pText = *s->Legend() + TString(" (#int L = ") + TString(l) + TString("/pb)");
512			TString pText = *s->Legend() + TString(" (") + TString(l) + TString("/pb)");
513			MDB(kGeneral,1)
514			printf(" Adding text \"%s\" at: (x,y) = (%6.2f,%6.2f)\n",pText.Data(),xText,yText);
515			// set the proper text color
516			text->SetTextColor(hStyle->Color());
517			// plot the text
518			text->SetTextAlign(12);
519			text->DrawLatex(xText,yText,pText.Data());
520			}
521			}
522
523			// Make sure the histogram frame is nice
524			box->SetFillStyle(0);
525			box->SetLineColor(kBlack);
526			box->DrawBox(0,0,100,100);
527
528			delete text;
529			delete box;
530
531			return;
532			}