BCH2D Class Reference

A class for handling 2D distributions. More...

#include <BCH2D.h>

List of all members.

Public Member Functions

Constructors and destructors



 BCH2D (TH2D *h)
 BCH2D ()
 ~BCH2D ()
Member functions (miscellaneous methods)



void CalculateIntegratedHistogram ()
void Draw (int options=0, bool drawmode=true)
TGraph * GetBandGraph (TH2D *h, double level1, double level2)
TGraph * GetBandGraph (double level1, double level2)
TGraph ** GetBandGraphs (TH2D *h, int &n)
double GetLevel (double p)
std::vector< double > GetLevelBoundary (TH2D *h, double level)
std::vector< double > GetLevelBoundary (double level)
TGraph * GetLowestBandGraph (TH2D *h)
TGraph * GetLowestBandGraph (TH2D *h, std::vector< int > nint)
std::vector< int > GetNIntervalsY (TH2D *h, int &nfoundmax)
void Print (const char *filename, int options=0, int ww=0, int wh=0)
Member functions (get)



TH2D * GetHistogram ()
void GetMode (double &mode)
Member functions (set)



void SetGlobalMode (double mode[2])
void SetHistogram (TH2D *hist)

Private Attributes

TH2D * fHistogram
TH1D * fIntegratedHistogram
double fMode [2]
int fModeFlag

Detailed Description

A class for handling 2D distributions.

Author:
Daniel Kollar
Kevin Kröninger
Version:
1.0
Date:
08.2008 This class contains a TH2D histogram and some additional functions. It is used for marginalized distributions.

Definition at line 33 of file BCH2D.h.


Constructor & Destructor Documentation

BCH2D::BCH2D (  ) 

The default constructor.

Definition at line 28 of file BCH2D.cxx.

00029 {
00030    fHistogram = 0;
00031    fIntegratedHistogram = 0;
00032 
00033    fModeFlag = 0;
00034 }

BCH2D::BCH2D ( TH2D *  h  ) 

The complete constructor.

Definition at line 38 of file BCH2D.cxx.

00039 {
00040    fHistogram = h;
00041    fIntegratedHistogram = 0;
00042 
00043    fModeFlag = 0;
00044 }

BCH2D::~BCH2D (  ) 

The default destructor.

Definition at line 48 of file BCH2D.cxx.

00049 {
00050    if (fHistogram) delete fHistogram;
00051    if (fIntegratedHistogram) delete fIntegratedHistogram;
00052 }


Member Function Documentation

void BCH2D::CalculateIntegratedHistogram (  ) 

Definition at line 345 of file BCH2D.cxx.

00346 {
00347    int nz = 100;
00348 
00349    double zmax = fHistogram -> GetMaximum();
00350    double dz   = zmax / double(nz);
00351 
00352    double nx = fHistogram -> GetNbinsX();
00353    double ny = fHistogram -> GetNbinsY();
00354 
00355    // create histogram
00356    if (fIntegratedHistogram)
00357       delete fIntegratedHistogram;
00358 
00359    fIntegratedHistogram = new TH1D(
00360          TString::Format("%s_int_prob_%d",fHistogram->GetName(),BCLog::GetHIndex()), "", nz, 0.0, zmax);
00361    fIntegratedHistogram -> SetXTitle("z");
00362    fIntegratedHistogram -> SetYTitle("Integrated probability");
00363    fIntegratedHistogram -> SetStats(kFALSE);
00364 
00365    // loop over histogram
00366    for (int ix = 1; ix <= nx; ix++)
00367       for (int iy = 1; iy <= ny; iy++)
00368       {
00369          int binmin = BCMath::Nint(fHistogram -> GetBinContent(ix, iy) / dz);
00370          for (int i = binmin; i <= nz; i++)
00371             fIntegratedHistogram -> SetBinContent(i,
00372                   fIntegratedHistogram -> GetBinContent(i) +
00373                   fHistogram -> GetBinContent(ix, iy));
00374       }
00375 }

void BCH2D::Draw ( int  options = 0,
bool  drawmode = true 
)

Draw 2-d distribution into the active canvas

Parameters:
options explanation to come
drawmode specify whether a marker should be drawn at the location of the mode

Definition at line 92 of file BCH2D.cxx.

00093 {
00094    // draw histogram
00095    fHistogram -> SetLineColor(kBlack);
00096 // fHistogram -> SetLineWidth(4);
00097 
00098    fHistogram->Scale(1./fHistogram->Integral("width"));
00099 
00100    double modex,modey;
00101 
00102    // set mode to display
00103    if(fModeFlag)
00104    {
00105       modex = fMode[0];
00106       modey = fMode[1];
00107    }
00108    else
00109    {
00110       int maximumbin = fHistogram -> GetMaximumBin();
00111 
00112       int binx = maximumbin % (fHistogram -> GetNbinsX() + 2);
00113       int biny = maximumbin / (fHistogram -> GetNbinsX() + 2);
00114 
00115       modex = fHistogram -> GetXaxis() -> GetBinCenter(binx);
00116       modey = fHistogram -> GetYaxis() -> GetBinCenter(biny);
00117    }
00118 
00119    // normalize histogram
00120    fHistogram->Scale(1./fHistogram->Integral("width"));
00121 
00122    // draw according to selected option
00123    if (options == 0)
00124       fHistogram -> Draw("CONT0");
00125    else if (options == 1)
00126    {
00127       fHistogram -> Draw("CONT3");
00128 
00129       // set contours
00130       this -> CalculateIntegratedHistogram();
00131 
00132       double levels[4];
00133       levels[0] = 0.;
00134       levels[1] = this -> GetLevel(1.0 - 0.6827);
00135       levels[2] = this -> GetLevel(1.0 - 0.9545);
00136       levels[3] = this -> GetLevel(1.0 - 0.9973);
00137 
00138       fHistogram -> SetContour(4, levels);
00139 
00140       // best fit value
00141       TMarker* marker = new TMarker(modex, modey, 24);
00142       marker -> Draw();
00143 
00144       TLegend* legend = new TLegend(0.65, 0.80, 0.95, 0.95);
00145       legend -> SetBorderSize(0);
00146       legend -> SetFillColor(kWhite);
00147       legend -> AddEntry(fHistogram, "68% prob. region", "L");
00148       legend -> AddEntry(marker, "Best fit", "P");
00149       legend -> Draw();
00150    }
00151    else if (options == 2)
00152    {
00153       fHistogram -> Draw("CONT3");
00154 
00155       // set contours
00156       this -> CalculateIntegratedHistogram();
00157 
00158       double levels[2];
00159       double level32 = this -> GetLevel(0.32);
00160       levels[0] = 0.;
00161       levels[1] = level32;
00162 
00163       fHistogram -> SetContour(2, levels);
00164 
00165       // best fit value
00166       TMarker* marker = new TMarker(modex, modey, 24);
00167       marker -> Draw();
00168 
00169       TLegend* legend = new TLegend(0.65, 0.80, 0.95, 0.95);
00170       legend -> SetBorderSize(0);
00171       legend -> SetFillColor(kWhite);
00172       legend -> AddEntry(fHistogram, "68% prob. region", "L");
00173       legend -> AddEntry(marker, "Best fit", "P");
00174       legend -> Draw();
00175 
00176    }
00177    else if (options == 3)
00178    {
00179       fHistogram -> Draw("CONT3");
00180 
00181       // set contours
00182       this -> CalculateIntegratedHistogram();
00183 
00184       double levels[2];
00185       double level10 = this -> GetLevel(0.10);
00186       levels[0] = 0.;
00187       levels[1] = level10;
00188 
00189       fHistogram -> SetContour(2, levels);
00190 
00191       TLegend* legend = new TLegend(0.65, 0.80, 0.95, 0.95);
00192       legend -> SetBorderSize(0);
00193       legend -> SetFillColor(kWhite);
00194       legend -> AddEntry(fHistogram, "90% prob. region", "L");
00195       legend -> Draw();
00196    }
00197    else if (options == 4)
00198    {
00199       fHistogram -> Draw("CONT3");
00200 
00201       // set contours
00202       this -> CalculateIntegratedHistogram();
00203 
00204       double levels[2];
00205       double level5 = this -> GetLevel(0.05);
00206       levels[0] = 0.;
00207       levels[1] = level5;
00208 
00209       fHistogram -> SetContour(2, levels);
00210 
00211       TLegend* legend = new TLegend(0.65, 0.80, 0.95, 0.95);
00212       legend -> SetBorderSize(0);
00213       legend -> SetFillColor(kWhite);
00214       legend -> AddEntry(fHistogram, "95% prob. region", "L");
00215       legend -> Draw();
00216    }
00217    else if (options == 5)
00218       fHistogram -> Draw("COL");
00219    else if (options == 52 || options == 521)
00220    {
00221       // create new empty histogram
00222       int nx = fHistogram -> GetNbinsX();
00223       int ny = fHistogram -> GetNbinsY();
00224       double xmin = fHistogram->GetXaxis()->GetXmin();
00225       double xmax = fHistogram->GetXaxis()->GetXmax();
00226       TH2D * h = new TH2D(
00227             TString::Format("htemp52_%d",BCLog::GetHIndex()),fHistogram->GetTitle(),
00228             nx,xmin,xmax,
00229             ny,fHistogram->GetYaxis()->GetXmin(),fHistogram->GetYaxis()->GetXmax());
00230       h->SetXTitle(fHistogram->GetXaxis()->GetTitle());
00231       h->SetYTitle(fHistogram->GetYaxis()->GetTitle());
00232 
00233       // copy contents of the main histogram
00234 //    double min = fHistogram->GetMinimum(0.);
00235       for(int i=1;i<=nx;i++)
00236          for(int j=1;j<=ny;j++)
00237          {
00238             double val = fHistogram -> GetBinContent(i,j);
00239             // if requested, change contents of bins to log scale
00240             if(options == 521)
00241             {
00242 //             if(val == 0.)
00243 //                val = log(min)-1.;
00244 //             else
00245                   val = log10(val);
00246             }
00247             h -> SetBinContent(i,j,val);
00248          }
00249 
00250       // draw
00251       h->SetStats(0);
00252       h->Draw("col");
00253 
00254       // draw contour
00255       fHistogram -> Draw("cont3 same");
00256       fHistogram -> SetLineWidth(2);
00257 
00258       // set contours
00259       this -> CalculateIntegratedHistogram();
00260 
00261       double levels[] = { this -> GetLevel(0.32) };
00262       fHistogram -> SetContour(1, levels);
00263 
00264       // best fit value
00265       if(drawmode)
00266       {
00267          TMarker * marker0 = new TMarker(modex, modey, 8);
00268          marker0 -> SetMarkerColor(0);
00269          marker0 -> SetMarkerSize(.7);
00270          marker0 -> Draw();
00271          TMarker * marker1 = new TMarker(modex, modey, 4);
00272          marker1 -> SetMarkerColor(1);
00273          marker1 -> SetMarkerSize(.7);
00274          marker1 -> Draw();
00275 //       TMarker * marker = new TMarker(modex, modey, 5);
00276 //       marker -> SetMarkerColor(0);
00277 //       marker -> Draw();
00278       }
00279    }
00280    else if (options == 53 || options == 531)
00281    {
00282       // create new empty histogram
00283 //    int nx = fHistogram -> GetNbinsX();
00284 //    int ny = fHistogram -> GetNbinsY();
00285 //    double xmin = fHistogram->GetXaxis()->GetXmin();
00286 //    double xmax = fHistogram->GetXaxis()->GetXmax();
00287 //    TH2D * h = new TH2D(
00288 //          TString::Format("htemp53_%d",BCLog::GetHIndex()),fHistogram->GetTitle(),
00289 //          nx,xmin,xmax,
00290 //          ny,fHistogram->GetYaxis()->GetXmin(),fHistogram->GetYaxis()->GetXmax());
00291 //    h->SetXTitle(fHistogram->GetXaxis()->GetTitle());
00292 //    h->SetYTitle(fHistogram->GetYaxis()->GetTitle());
00293 
00294       // copy contents of the main histogram
00295 //    double min = fHistogram->GetMinimum(0.);
00296 /*    for(int i=1;i<=nx;i++)
00297          for(int j=1;j<=ny;j++)
00298          {
00299             double val = fHistogram -> GetBinContent(i,j);
00300             // if requested, change contents of bins to log scale
00301             if(options == 531)
00302             {
00303 //             if(val == 0.)
00304 //                val = log(min)-1.;
00305 //             else
00306                   val = log10(val);
00307             }
00308             h -> SetBinContent(i,j,val);
00309          }
00310 
00311       // draw
00312       h->SetStats(0);
00313       h->Draw("colz");
00314 */
00315       gPad -> SetLogz();
00316       fHistogram -> Draw("colz");
00317 
00318       // draw contour
00319 //    fHistogram -> Draw("cont3 same");
00320 //    fHistogram -> SetLineWidth(2);
00321 
00322       // set contours
00323 //    this -> CalculateIntegratedHistogram();
00324 
00325 //    double levels[] = { this -> GetLevel(0.32) };
00326 //    fHistogram -> SetContour(1, levels);
00327 
00328       // best fit value
00329       TMarker * marker0 = new TMarker(modex, modey, 8);
00330       marker0 -> SetMarkerColor(0);
00331       marker0 -> SetMarkerSize(.7);
00332       marker0 -> Draw();
00333       TMarker * marker1 = new TMarker(modex, modey, 4);
00334       marker1 -> SetMarkerColor(1);
00335       marker1 -> SetMarkerSize(.7);
00336       marker1 -> Draw();
00337 //    TMarker * marker = new TMarker(modex, modey, 5);
00338 //    marker -> SetMarkerColor(0);
00339 //    marker -> Draw();
00340    }
00341 }

TGraph * BCH2D::GetBandGraph ( TH2D *  h,
double  level1,
double  level2 
)

Definition at line 562 of file BCH2D.cxx.

00563 {
00564    // define new graph
00565    int nx = h -> GetNbinsX() - 1;
00566 
00567    TGraph * g = new TGraph(2*nx);
00568 // g -> SetFillStyle(1001);
00569 // g -> SetFillColor(kYellow);
00570 
00571    // get error bands
00572    std::vector <double> ymin = GetLevelBoundary(h,l1);
00573    std::vector <double> ymax = GetLevelBoundary(h,l2);
00574 
00575    for (int i = 0; i < nx; i++)
00576    {
00577       g -> SetPoint(i, h -> GetXaxis() -> GetBinCenter(i+1), ymin[i]);
00578       g -> SetPoint(nx+i, h -> GetXaxis() -> GetBinCenter(nx-i), ymax[nx-i-1]);
00579    }
00580 
00581    return g;
00582 }

TGraph * BCH2D::GetBandGraph ( double  level1,
double  level2 
)

Definition at line 555 of file BCH2D.cxx.

00556 {
00557    return GetBandGraph(fHistogram , l1, l2);
00558 }

TGraph ** BCH2D::GetBandGraphs ( TH2D *  h,
int &  n 
)

Definition at line 392 of file BCH2D.cxx.

00393 {
00394    n=0;
00395 
00396    int nbands=0;
00397    TH2D * hcopy = (TH2D*)h->Clone(TString::Format("%s_copy_%d",h->GetName(),BCLog::GetHIndex()));
00398 
00399    std::vector <int> nint=GetNIntervalsY(hcopy,nbands);
00400 
00401    if(nbands>2)
00402    {
00403       BCLog::Out(BCLog::warning,BCLog::warning,
00404          Form("BCH2D::GetBandGraphs : Detected %d bands. Maximum allowed is 2 (sorry).",nbands));
00405       return 0;
00406    }
00407    else if(nbands==0)
00408    {
00409       BCLog::Out(BCLog::warning,BCLog::warning,"BCH2D::GetBandGraphs : No bands detected.");
00410       return 0;
00411    }
00412 
00413    TGraph ** gxx = new TGraph*[nbands];
00414 
00415    TH2D * h0 = (TH2D*)h->Clone();
00416 
00417    if (nbands>0)
00418       gxx[0] = GetLowestBandGraph(h0,nint);
00419 
00420    if (nbands==2)
00421    {
00422       gxx[1] = GetLowestBandGraph(h0);
00423       n=2;
00424    }
00425    else
00426       n=1;
00427 
00428    return gxx;
00429 }

TH2D* BCH2D::GetHistogram (  )  [inline]
Returns:
The 2D histogram.

Definition at line 59 of file BCH2D.h.

00060          { return fHistogram; };

double BCH2D::GetLevel ( double  p  ) 

Definition at line 379 of file BCH2D.cxx.

00380 {
00381    double quantiles[1];
00382    double probsum[1];
00383    probsum[0] = p;
00384 
00385    fIntegratedHistogram -> GetQuantiles( 1, quantiles, probsum);
00386 
00387    return quantiles[0];
00388 }

std::vector< double > BCH2D::GetLevelBoundary ( TH2D *  h,
double  level 
)

Definition at line 528 of file BCH2D.cxx.

00529 {
00530    std::vector <double> b;
00531 
00532    int nx = h -> GetNbinsX();
00533 
00534    b.assign(nx - 1, 0.0);
00535 
00536    // loop over x and y bins.
00537    for (int ix = 1; ix < nx; ix++)
00538    {
00539       TH1D * h1 = h -> ProjectionY(TString::Format("temphist_%d",BCLog::GetHIndex()), ix, ix);
00540 
00541       int nprobSum = 1;
00542       double q[1];
00543       double probSum[] = { level };
00544 
00545       h1 -> GetQuantiles(nprobSum, q, probSum);
00546 
00547       b[ix-1] = q[0];
00548    }
00549 
00550    return b;
00551 }

std::vector< double > BCH2D::GetLevelBoundary ( double  level  ) 

Definition at line 521 of file BCH2D.cxx.

00522 {
00523    return GetLevelBoundary(fHistogram, level);
00524 }

TGraph * BCH2D::GetLowestBandGraph ( TH2D *  h  ) 

Definition at line 469 of file BCH2D.cxx.

00470 {
00471    int n;
00472    return GetLowestBandGraph(h,GetNIntervalsY(h,n));
00473 }

TGraph * BCH2D::GetLowestBandGraph ( TH2D *  h,
std::vector< int >  nint 
)

Definition at line 477 of file BCH2D.cxx.

00478 {
00479    int nx = h -> GetNbinsX() - 1;
00480    int ny = h -> GetNbinsY();
00481 
00482    TGraph * g = new TGraph(2*nx);
00483 
00484    for (int ix=1; ix<=nx; ix++)
00485    {
00486       // get x for the bin
00487       double x;
00488       if(ix==1)
00489          x = h->GetXaxis()->GetBinLowEdge(1);
00490       else if(ix==nx)
00491          x = h->GetXaxis()->GetBinLowEdge(nx+1);
00492       else
00493          x = h->GetXaxis()->GetBinCenter(ix);
00494 
00495       for(int iy=1; iy<=ny; iy++)
00496          if(h->GetBinContent(ix,iy)>0.)
00497          {
00498             // get low edge of the first not empty bin in y
00499             g->SetPoint(ix-1, x, h->GetYaxis()->GetBinLowEdge(iy));
00500 
00501             // delete content of all subsequent not empty bins
00502             if(nint[ix-1]==2)
00503                h->SetBinContent(ix,iy,0.);
00504 
00505             while(h->GetBinContent(ix,++iy)>0.)
00506                if(nint[ix-1]==2)
00507                   h->SetBinContent(ix,iy,0.);
00508 
00509             // get low edge of the first empty bin in y
00510             g->SetPoint(2*nx-ix, x, h->GetYaxis()->GetBinLowEdge(iy));
00511 
00512             break;
00513          }
00514    }
00515 
00516    return g;
00517 }

void BCH2D::GetMode ( double &  mode  ) 
Returns:
The mean of the distribution.
The mode of the distribution

Definition at line 56 of file BCH2D.cxx.

00057 {
00058    // what is this?
00059 
00060 // double mode[2];
00061 // int binx, biny, binz;
00062 
00063 // fHistogram -> GetMaximumBin(binx, biny, binz);
00064 // mode = fHistogram -> GetBinCenter();
00065 
00066 // return mode;
00067 }

std::vector< int > BCH2D::GetNIntervalsY ( TH2D *  h,
int &  nfoundmax 
)

Returns the number of intervals as a function of x

Parameters:
h The histogram.
nfoundmax The maximum number of intervals.
Returns:
A vector containing the number of intervals for all bins in x.

Definition at line 433 of file BCH2D.cxx.

00434 {
00435    std::vector <int> nint;
00436 
00437    int nx = h -> GetNbinsX();
00438    int ny = h -> GetNbinsY();
00439 
00440    nfoundmax=0;
00441 
00442    // loop over histogram bins in x
00443    for (int ix=1; ix<=nx; ix++)
00444    {
00445       int nfound=0;
00446 
00447       // loop over histogram bins in y
00448       // count nonzero intervals in y
00449       for (int iy=1; iy<=ny; iy++)
00450          if(h->GetBinContent(ix,iy)>0.)
00451          {
00452             while(h->GetBinContent(ix,++iy)>0.)
00453                ;
00454             nfound++;
00455          }
00456 
00457       // store maximum number of nonzero intervals for the histogram
00458       if(nfound>nfoundmax)
00459          nfoundmax=nfound;
00460 
00461       nint.push_back(nfound);
00462    }
00463 
00464    return nint;
00465 }

void BCH2D::Print ( const char *  filename,
int  options = 0,
int  ww = 0,
int  wh = 0 
)

Print 2-d histogram to file

Parameters:
filename The filename
ww canvas size in pixels along X
ww canvas size in pixels along Y If ww and wh are set to 0, default ROOT canvas size is used. For explanation of the parameter options see the Draw() method.

Definition at line 71 of file BCH2D.cxx.

00072 {
00073    // create temporary canvas
00074    TCanvas * c;
00075    if(ww >0 && wh > 0)
00076       c = new TCanvas("c","c",ww,wh);
00077    else
00078       c = new TCanvas("c","c",700,700);
00079    c -> cd();
00080 
00081    // draw histogram
00082    this->Draw(options);
00083 
00084    gPad->RedrawAxis();
00085 
00086    // print to file
00087    c -> Print(filename);
00088 }

void BCH2D::SetGlobalMode ( double  mode[2]  )  [inline]

Set global mode.

Parameters:
The global mode.

Definition at line 82 of file BCH2D.h.

00083          { fMode[0] = mode[0]; fMode[1] = mode[1]; fModeFlag =1; };

void BCH2D::SetHistogram ( TH2D *  hist  )  [inline]

Set the 2D histogram.

Definition at line 76 of file BCH2D.h.

00077          { fHistogram = hist; };


Member Data Documentation

TH2D* BCH2D::fHistogram [private]

The 2D histogram

Definition at line 142 of file BCH2D.h.

TH1D* BCH2D::fIntegratedHistogram [private]

The integrated 2D histogram

Definition at line 146 of file BCH2D.h.

double BCH2D::fMode[2] [private]

Global mode

Definition at line 150 of file BCH2D.h.

int BCH2D::fModeFlag [private]

"Is there a global mode?" flag

Definition at line 154 of file BCH2D.h.


The documentation for this class was generated from the following files:

Generated on Tue Oct 6 09:48:21 2009 for Bayesian Analysis Toolkit by  doxygen 1.6.1