CorrelatedSpaceTH2.cpp

Go to the documentation of this file.

/***************************************************************************
 *   Copyright (C) 2004-2006 by Olivier Stezowski & Christian Finck        *
 *   stezow(AT)ipnl.in2p3.fr, cfinck(AT)ires.in2p3.fr                      *
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 *   This program is distributed in the hope that it will be useful,       *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
 *   GNU General Public License for more details.                          *
 *                                                                         *
 *   You should have received a copy of the GNU General Public License     *
 *   along with this program; if not, write to the                         *
 *   Free Software Foundation, Inc.,                                       *
 *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
 ***************************************************************************/

#ifndef Gw_CorrelatedSpaceTH2
#include "CorrelatedSpaceTH2.h"
#endif

#ifndef Gw_CSFilter
#include "CSFilter.h"
#endif

#ifndef ROOT_TBrowser
#include "TBrowser.h"
#endif

#ifndef ROOT_TFile
#include "TFile.h"
#endif

#ifndef ROOT_TH2
#include "TH2.h"
#endif

#ifndef ROOT_TROOT
#include "TROOT.h"
#endif

#ifndef ROOT_TRandom
#include "TRandom.h"
#endif

#ifndef ROOT_TString
#include "TString.h"
#endif

#ifndef ROOT_TTree
#include "TTree.h"
#endif

#ifndef ROOT_TVirtualPad
#include "TVirtualPad.h"
#endif


using namespace Gw;

//
ClassImp(CorrelatedSpaceTH2)

//__________________________________________________________
CorrelatedSpaceTH2::CorrelatedSpaceTH2(const char* name, const char* title, Option_t* type) 
  : ASymCorrelatedSpace(name, title),
        fHisto(0x0),
    fAdopted(false),
    fLog("CorrelatedSpaceTH2")
{ 
  // default constructor

  SetDimension(2);
  TString sType(type);
  sType.ToLower();
  if (sType.Contains("float"))
        fType = 0;
  else if (sType.Contains("short"))
        fType = 1;
  else if (sType.Contains("int"))
        fType = 2;
  else if (sType.Contains("double"))
        fType = 3;
  else {
        fLog.SetProcessMethod("CorrelatedSpaceTH2(Option_t* )");
        fLog << warning << "Unkwon type, set to float" << dolog;
        fType = 0;
  }
}

//__________________________________________________________
CorrelatedSpaceTH2::~CorrelatedSpaceTH2()
{
  // default destructor
   if (!fAdopted)
          delete fHisto;
}

//__________________________________________________________
const Char_t* CorrelatedSpaceTH2::GetHistoName() const
{
  return Form("%s_H2", GetName()); 
}

//__________________________________________________________
const Char_t* CorrelatedSpaceTH2::GetHistoTitle() const
{
  return Form("%s_H2", GetTitle()); 
}

//__________________________________________________________
Double_t CorrelatedSpaceTH2::GetContent(Int_t i, Int_t j)
{
  
  if (fHisto == 0x0) {
        fLog.SetProcessMethod("GetContent(Int_t, Int_t )");
        fLog << error << "Histo TH2 not definied" << dolog;
        return 0.;
  }
  
  return (Double_t)fHisto->GetBinContent(fHisto->GetBin(i,j));
}

//__________________________________________________________
Double_t CorrelatedSpaceTH2::GetCellContent(Int_t iCell)
{
  if (fHisto == 0x0) {
        fLog.SetProcessMethod("GetCellContent(Int_t )");
        fLog << error << "Histo TH2 not definied" << dolog;
        return 0.;
  }
  
  Int_t dimY = fHisto->GetNbinsY();
  
  Int_t j = iCell / dimY;
  Int_t i = iCell % dimY;

  return GetContent(i, j);
}

//__________________________________________________________
Int_t CorrelatedSpaceTH2::GetNbin(Int_t axis) const
{
  if (fHisto == 0x0) {
        fLog.SetProcessMethod("GetNbin(Int_t )");
        fLog << error << "Histo TH2 not definied" << dolog;
        return -1;
  }
  
  switch (axis) {
        case 0:
          return fHisto->GetNbinsX();
        case 1:
          return fHisto->GetNbinsY();
        default:
          return -1;
  }
}

//__________________________________________________________
void CorrelatedSpaceTH2::SetCellContent(Int_t iCell, Double_t weight)
{
  AssertHisto();
  Int_t dimY = GetHisto("y")->GetNbinsX();

  Int_t j = iCell / (dimY);
  Int_t i = iCell % (dimY);
  
  fHisto->Fill(i, j, weight);
  
  GetHistoTot("x")->Fill(j, weight);
  GetHistoTot("y")->Fill(i, weight);
}

//__________________________________________________________
void CorrelatedSpaceTH2::AssertHisto()
{
  if (fHisto == 0x0) {
        switch (fType) {
          case 0:
                fHisto = new TH2F();
                break;
          case 1:
                fHisto = new TH2S();
                break;
          case 2:
                fHisto = new TH2I();
                break;
          case 3:
                fHisto = new TH2D();
                break;
        }
        fHisto->SetName(GetHistoName());
        fHisto->SetTitle(GetHistoTitle());
        Int_t dimX   = GetHisto("x")->GetNbinsX();
        Int_t dimY   = GetHisto("y")->GetNbinsX();
        Float_t minX = GetHisto("x")->GetXaxis()->GetXmin();
        Float_t maxX = GetHisto("x")->GetXaxis()->GetXmax();
        Float_t minY = GetHisto("y")->GetXaxis()->GetXmin();
        Float_t maxY = GetHisto("y")->GetXaxis()->GetXmax();
        fHisto->SetBins(dimX, minX, maxX, dimY, minY, maxY);
  }
}

//__________________________________________________________
void CorrelatedSpaceTH2::Fill(const Double_t* Xn, const Double_t weight)
{
  
  AssertHisto();
  
  Int_t dimX = GetHisto("x")->GetNbinsX();
  Int_t dimY = GetHisto("y")->GetNbinsX();
 
   
  fHisto->SetBins(dimX, 0, dimX, dimY, 0, dimY);
  SetNofCells(dimX*dimY);
  
  for (Int_t i = 0; i < dimX; ++i) {
        for (Int_t j = 0; j < dimY; ++j) {
          fHisto->Fill(i,j, Xn[i*dimX + j]*weight);
        }       
  }
  
  TH1D* h1x = fHisto->ProjectionX();
  TH1D* h1y = fHisto->ProjectionY();

  for (Int_t i= 0; i < dimX; ++i) 
        GetHistoTot("x")->Fill(i, h1x->GetBinContent(i));
  
  for (Int_t j = 0; j < dimY; ++j)  
        GetHistoTot("y")->Fill(j, h1y->GetBinContent(j));
}

//__________________________________________________________
void CorrelatedSpaceTH2::FillFromH2(TH2* h2)
{
    
    fLog.SetProcessMethod("FillFromH2(TH2* )");
    
    if (gFile == 0x0) {
        fLog << error << "No file exists to write" << dolog;
        return;
    }
    
    Int_t dimX = h2->GetNbinsX();
    Int_t dimY = h2->GetNbinsY();
    Float_t minX = h2->GetXaxis()->GetXmin(); // modif des 4 lignes
    Float_t maxX = h2->GetXaxis()->GetXmax();
    Float_t minY = h2->GetXaxis()->GetXmin();
    Float_t maxY = h2->GetXaxis()->GetXmax();
    
    SetBins(dimX, minX, maxX, 0);
    SetBins(dimY, minY, maxY, 1);
    
    SetNofCells(dimX*dimY);
    
    fHisto = (TH2*)h2->Clone();
    TH1D* h1x = h2->ProjectionX();
    TH1D* h1y = h2->ProjectionY();
    
    for (Int_t i= 0; i < dimX; ++i)
        GetHistoTot("x")->SetBinContent(i, h1x->GetBinContent(i));
        
    for (Int_t j = 0; j < dimY; ++j)
        GetHistoTot("y")->SetBinContent(j, h1y->GetBinContent(j));
    
    fHisto->SetName(GetHistoName());
    fHisto->SetTitle(GetHistoTitle());
    
    fLog << dolog;
}


//__________________________________________________________
void CorrelatedSpaceTH2::FillRandom(Int_t dimX, Int_t dimY)
{
  fLog.SetProcessMethod("FillRandom(Int_t )");
  
  AssertHisto();
  
  SetBins(dimX, 0, dimX, 0);
  SetBins(dimY, 0, dimY, 1);
  
  fHisto->SetBins(dimX, 0, dimX, dimY, 0, dimY);

  SetNofCells(dimX*dimY);

  if (gFile == 0x0) {
        fLog << error << "No file exists to write" << dolog;
        return;
  }
  
  Double_t counts = 0;
  for (Int_t i = 0; i < dimX; ++i) {
        printf("%d\r", i);
        for (Int_t j = 0; j < dimY; ++j) {
          counts = gRandom->Uniform(0, 10);
          fHisto->Fill(i,j, counts);
        }       
  }
  
  TH1D* h1x = fHisto->ProjectionX();
  TH1D* h1y = fHisto->ProjectionY();
  
  for (Int_t i= 0; i < dimX; ++i) 
        GetHistoTot("x")->Fill(i, h1x->GetBinContent(i));
  
  for (Int_t j = 0; j < dimY; ++j)  
        GetHistoTot("y")->Fill(j, h1y->GetBinContent(j));

  fLog << dolog; 
}

//__________________________________________________________
void CorrelatedSpaceTH2::Adopt(TH2* h2)
{
   if (fHisto) delete fHisto;
   fHisto = h2;
   
   Int_t dimX = h2->GetNbinsX();
   Float_t minX = h2->GetXaxis()->GetXmin();
   Float_t maxX = h2->GetXaxis()->GetXmax();
   SetBins(dimX, minX, maxX, 0);
   
   Int_t dimY = h2->GetNbinsY();
   Float_t minY = h2->GetYaxis()->GetXmin();
   Float_t maxY = h2->GetYaxis()->GetXmax();
   SetBins(dimY, minY, maxY, 1);
   
   SetNofCells(dimX*dimY);
   
   TH1D* h1x = fHisto->ProjectionX();
   TH1D* h1y = fHisto->ProjectionY();
   h1x->SetTitle("");
   h1y->SetTitle("");
   
   for (Int_t i= 0; i < dimX; ++i) 
          GetHistoTot("x")->SetBinContent(i, h1x->GetBinContent(i));
   
   for (Int_t j = 0; j < dimY; ++j)  
          GetHistoTot("y")->SetBinContent(j, h1y->GetBinContent(j));
   
   fAdopted = true;
}

//__________________________________________________________
void CorrelatedSpaceTH2::Project(CSFilter* filter, Option_t* axis)
{  
  fLog.SetProcessMethod("Project(CSFilter*, Option_t* )");
  
  if (fStatus == kUndefined) {
        fLog << error << "TH2 not defined" << dolog;
        return;
  }
  
  if (filter == 0x0) {
        fLog << error << "No filter definied yet" << dolog;
        return;
  }
  
  filter->SetCondition(1);
  filter->SetDimension(fDimension-1);
  
  Int_t idx = GetAxisNumber(axis);
  
  if (idx > 1) {
        fLog << warning << "Axis number out of range" << dolog;
        return;
  }
  
  TH1D* h1     = 0x0;
  TH1D* hProj1 = 0x0;  
  TH1D* hProj2 = 0x0;

  if ( filter->IsForBkg() )
        h1 =  static_cast<TH1D*> ( fHistoBkgList->At(idx) );
  else 
        h1 = static_cast<TH1D*> ( fHistoList->At(idx) );  
  
  h1->Reset();
  
  Int_t dim1 = 0;
  Int_t dim2 = 0;
  Double_t counts = 0.;
  
  if (idx == 0) {
        dim1 = GetHistoTot("x")->GetNbinsX();
        dim2 = GetHistoTot("y")->GetNbinsX();
        hProj1 = GetHistoTot("x");
        hProj2 = GetHistoTot("y");
  } else { 
        dim2 = GetHistoTot("x")->GetNbinsX();
        dim1 = GetHistoTot("y")->GetNbinsX();
        hProj2 = GetHistoTot("x");
        hProj1 = GetHistoTot("y");
  }

  Int_t    nofGates;
  Int_t    dimProj = fDimension - 1;
  Short_t  inGate[dimProj];
  Double_t energies[dimProj];
  
  for (Int_t i= 0; i < dim1; ++i) {
        
        energies[0] = hProj1->GetBinCenter(i);

        nofGates = filter->IsInside(energies, dimProj, inGate);

        if (nofGates > 1) {
          fLog << error << Form("Too much good energies : %d", nofGates) << dolog;
          return;
        }
        
        if (inGate[0] == 1) {
          for (Int_t j = 0; j < dim2; ++j) { 
                if (idx == 0)
                  counts = (Double_t)fHisto->GetBinContent(fHisto->GetBin(i,j));
                else
                  counts = (Double_t)fHisto->GetBinContent(fHisto->GetBin(j,i));
                h1->Fill(hProj2->GetBinCenter(j), counts);        
          }
        }
  }
  fLog << dolog; 
}
  
//_____________________________________________________________________________
void CorrelatedSpaceTH2::Project(Option_t* axis, Option_t* gateName) 
{ 
  CorrelatedSpace::Project(axis, gateName); 
}

//______________________________________________________________________________
void CorrelatedSpaceTH2::Browse(TBrowser* b)
{
   fHisto->Browse(b);
}