THistoManager.C

Go to the documentation of this file.

#ifndef __THISTOMANAGER_H
#include "Marana/THistoManager.h" 
#endif


void THistoManager::ReadParameters(string parameters_file,string Tag,const Short_t SysNum) 
{  

 string ClassName="";
 
 list<string> PartTags;
 ReadParameters(parameters_file,Tag,ClassName,PartTags,SysNum);

}


void THistoManager::ReadParameters(string parameters_file,string Tag,string ClassName,
                                        list<string>& PartTags,const Short_t SysNum) 
{  
  string item;
//... open the Parameters file
  
  try {

   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 40 "<< endl;              
    string mess=error.str();      
     
    throw mess; 
    return;
   }



   ifstream inParams(parameters_file.c_str());
   
   if(!inParams)
     throw string("THistoManager::ReadParameters = Error in opening histogram parameter file! - stop");
   
   string end="end_";
   string start=Tag+":";
   
   bool stop=false;
   string token,type,readtag,dummy;
   
//... stream buffer with histogram definitions
   istringstream* parstrstr=0;
   string parstr;
   
    while (getline(inParams, item)) {
      if (item.find(start) == string::npos) { continue; }
      cout << endl << "Start reading loop " << Tag << ": " << endl;
      string bufx(item.substr(item.find(":") + 1));
      string::size_type more = 0;
      parstr = "";
      do {
        //... eat whitespace
        while ( bufx[0] == ' ' ) {
          bufx = bufx.substr(1,bufx.size());
        }
        //cout << bufx << endl;

        // check for multiple entries separated by ,      
        more = bufx.find(",");
        //cout << more << endl;
        if (more == string::npos) {
          token = bufx;
        } else {
          token = bufx.substr(0,more);          
          bufx = bufx.substr(more+1,bufx.size());
        }

        //cout << "token: " <<token << endl;
        //... optional histogram id offset in parentheses  
        Int_t hoffset = 0;
        string::size_type s1 = token.find("(")+1;
        string::size_type s2 = token.find(")")-s1;
        string oc = token.substr(s1,s2);  
        hoffset = atoi(oc.data());

        //cout << hoffset << endl;
        if (hoffset) {
          type = token.substr(0,s1 - 1).c_str();
        } else {
          type = token;
          while ( bufx[type.size()-1] == ' ' ) {//... eat whitespace
            type = type.substr(0,type.size()-1);
          }
         }       
        if(ClassName!="") type=type+ClassName;

        cout << endl 
          << "________________________________________________________" << endl
          << "type: '" << type <<"' histogram id offset: " << hoffset   << endl
          << "________________________________________________________" << endl
          << endl;

        //... on first iteration read histogram definitions into a buffer
        if ( parstr == "" ) {        
          parstr = "";
          while (!stop && getline(inParams,item)) {
            parstr += item;            
            parstr += '\n';
            if (item == end + Tag) { stop = true; }                       
          }
        }

        //cout << parstr << endl;
        parstrstr = new istringstream(parstr.data());
        ParseHistos(*parstrstr,type,hoffset,Tag,PartTags,SysNum);
        delete parstrstr;

        //cout << "token: " <<token << endl;
        //if (more!=string::npos) cout  << "_______________________" << endl;

      } while (more!=string::npos);

      stop = false;      
    } //... while



   cout << endl <<"--> number of h1f histos "<< fHisto1D[SysNum]->GetEntries()<<endl;
   cout << "--> number of h2f histos "<< fHisto2D[SysNum]->GetEntries()<<endl;

//... close the Parameters file
   inParams.close(); 

  } //... try 
  catch(string error){
     cout << error <<endl;;
     exit(0); 
  } catch(...){
     throw string("THistoManager::Parameters = Undefined error in reading parameters_file! - stop");
     exit(0); 
  }     
     
}



bool THistoManager::ParseHistos(
  istringstream& inParams, string&type, Int_t hoffset,
  string Tag, list<string>& PartTags,const Short_t SysNum)
{   
//     Read and parse histo definitions
    
    bool   stop=false;
    string item;
    string end="end_";

    while (!stop && inParams >> item) {
//...read and create h1f histos ;
        if(item=="h1f") {
                ReadParamsh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="h1f_sys"){
                ReadParamsh1f(inParams,type,hoffset,SysNum,kTRUE);
        }
        else if(item=="h2f"){
                ReadParamsh2f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="h2f_sys"){
                ReadParamsh2f(inParams,type,hoffset,SysNum,kTRUE);
        }
//...read and create h1f generic histos ;
        else if(item=="h1fgeneric") {
                ReadParamsh1fgeneric(inParams,type,hoffset,PartTags,SysNum,kFALSE);
        }
        else if(item=="h1fgeneric_sys") {
                ReadParamsh1fgeneric(inParams,type,hoffset,PartTags,SysNum,kTRUE);
        }

//... for bin histograms with fix bin width
        else if(item=="2dbin"){
                ReadParams2dbin(inParams,type,hoffset,SysNum);
        }
        else if(item=="3dbin"){
                ReadParams3dbin(inParams,type,hoffset,SysNum);
        }

//... for bin histograms with fix bin width
        else if(item=="1binfh1f"){
                ReadParams1binfh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="1binfh1f_sys"){
                ReadParams1binfh1f(inParams,type,hoffset,SysNum,kTRUE);
        }
//... for bin histogram with variable bin width
        else if(item=="1binvbh1f"){
                ReadParams1binVarSizeh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="1binvbh1f_sys"){
                ReadParams1binVarSizeh1f(inParams,type,hoffset,SysNum,kTRUE);
        }
//... for bin histograms with variable bin width
        else if(item=="1binvh1f"){
                ReadParams1binvh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="1binvh1f_sys"){
                ReadParams1binvh1f(inParams,type,hoffset,SysNum,kTRUE);
        }
//... for 2bin histograms with fix bin width
        else if(item=="2binfh1f"){
                ReadParams2binfh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="2binfh1f_sys"){
                ReadParams2binfh1f(inParams,type,hoffset,SysNum,kTRUE);
        }
//... for 2bin histograms with variable bin width
        else if(item=="2binvh1f"){
                ReadParams2binvh1f(inParams,type,hoffset,SysNum,kFALSE);
        }
        else if(item=="2binvh1f_sys"){
                ReadParams2binvh1f(inParams,type,hoffset,SysNum,kTRUE);
        } else if (item == end + Tag) {
            stop = true;
          }                     //if
          
        }                       //while 
  return stop;
}



 
THistoManager::THistoManager(){ 

    for(Int_t i=0;i<40;i++) {
        fHisto1D[i] = new TClonesArray("TH1F",1);
        fHisto2D[i] = new TClonesArray("TH2F",1);
        fHisto3D[i] = new TClonesArray("TH3F",1);
    }
   
}; 
THistoManager::THistoManager(string Parameters_file,string Tag){ 

    for(Int_t i=0;i<40;i++) {
        fHisto1D[i] = new TClonesArray("TH1F",1);
        fHisto2D[i] = new TClonesArray("TH2F",1);
        fHisto3D[i] = new TClonesArray("TH3F",1);
    }
    ReadParameters(Parameters_file,Tag);
};  
 
THistoManager::~THistoManager(){ 

//---- don't know why, but seems to block the jobs at the end (for 2D histos only ?)
//---- therefore very brutal solution: remove it ! crack !!
//  for(Int_t i=0;i<40;i++) {
// //...delete array if 1D histo   
//      if(fHisto1D[i]){
//        fHisto1D[i]->Delete();
//        delete fHisto1D[i];
//      }
//    
// //...delete array if 2D histo 
//      if(fHisto2D[i]){
//        fHisto2D[i]->Delete();
//        delete fHisto2D[i];
//      }  
//   }

};

void THistoManager::FillHisto(const int& num, const double& value, const double& wait,const Short_t SysNum){ 

//...fill histogram
   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 40 "<< endl;              
    string mess=error.str();      
    throw mess; 
    return;
   }
  
   
  if(num>fHisto1D[SysNum]->GetLast()) return;
   TH1F* histo= (TH1F*)fHisto1D[SysNum]->At(num);
   
   if(histo != NULL)histo->Fill(value,wait);
   else {
    if(SysNum==0) {
        ostringstream error;
        error << "THistoManager::FillHisto = Error 1D histogram num " 
              << num << " doesn't exist "<< endl;              
        string mess=error.str();      
        throw mess; 
    }
   }
}; 

void THistoManager::FillHisto(const int& num, const double& value_x, const double& value_y, const double& wait,const Short_t SysNum){ 

//...fill histogram h2f
   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 40 "<< endl;              
    string mess=error.str();      
     
    throw mess; 
    return;
   }

   if(num>fHisto2D[SysNum]->GetLast()) return;
   TH2F* histo= (TH2F*)fHisto2D[SysNum]->At(num);
   
   if(histo != NULL)histo->Fill(value_x,value_y,wait);
   else {
    if(SysNum==0) {
        ostringstream error;
        error << "THistoManager::FillHisto = Error 2D histogram num " 
              << num << " doesn't exist "<< endl;              
        string mess=error.str();      
        throw mess; 
    }
  
   }
}; 

void THistoManager::FillHisto(const int& num, const double& value_x, const double& value_y, const double& value_z, const double& wait,const Short_t SysNum)
{ 
//...fill histogram h3f
   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 40 "<< endl;              
    string mess=error.str();      
     
    throw mess; 
    return;
   }

   if(num>fHisto3D[SysNum]->GetLast()) return;
   TH3F* histo= (TH3F*)fHisto3D[SysNum]->At(num);
   
   if(histo != NULL)histo->Fill(value_x,value_y,value_z,wait);
   else {
    if(SysNum==0) {
        ostringstream error;
        error << "THistoManager::FillHisto = Error 3D histogram num " 
              << num << " doesn't exist "<< endl;              
        string mess=error.str();      
        throw mess; 
    }
  
   }
}; 


void THistoManager::Fill1BinHisto(Int_t hnum,Double_t binvalue,Double_t value,Double_t wait,const Short_t SysNum)
{ 
//     Fill 1dim. bin histograms 
   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 40 "<< endl;              
    string mess=error.str();      
    throw mess; 
    return;
   }

 if(hnum>fHisto1D[SysNum]->GetLast()) return;
 TH1F* basehisto= (TH1F*)fHisto1D[SysNum]->At(hnum);
   
 if(basehisto) {
        Double_t bmin,bmax;
        for(Int_t i=1;i<=basehisto->GetBinContent(1);i++) {
                bmin=basehisto->GetBinContent(1+i);
                bmax=basehisto->GetBinContent(2+i);
//              printf("Min = %f Max = %f \n",bmin,bmax);
                if(binvalue>bmin&&binvalue<bmax)
                        ((TH1F*)fHisto1D[SysNum]->At(hnum+i))->Fill(value,wait); 
        }
 } else {
    if(SysNum==0) {
       ostringstream error;
       error << "THistoManager::Fill1BinHisto = Error bin histograms hnum " 
             << hnum << " doesn't exist "<< endl;              
       string mess=error.str();      
       throw mess; 
    }
 }
} 






void THistoManager::Fill2BinHisto(Int_t hnum,Double_t binvalue,Double_t binvalue2,Double_t value,Double_t wait,const Short_t SysNum)
{ 
//     Fill 2dim. bin histograms 
   if(SysNum>=40) {
     ostringstream error;
      error << "THistoManager::FillHisto = Error Systematic num " 
          << SysNum << " >= 20 "<< endl;              
    string mess=error.str();      
     
    throw mess; 
    return;
   }

 if(hnum>fHisto1D[SysNum]->GetLast()) return;
 TH1F* basehisto= (TH1F*)fHisto1D[SysNum]->At(hnum);
   
 Int_t numhisto=hnum;  
 if(basehisto) {
        Double_t bmin,bmax,bmin2,bmax2;
//      cout<<basehisto->GetBinContent(1)<<"*"<<basehisto->GetBinContent(2)<<endl;
        for(Int_t i=1;i<=basehisto->GetBinContent(1);i++) {
           for(Int_t j=1;j<=basehisto->GetBinContent(2);j++) {
                bmin=basehisto->GetBinContent(2+i);
                bmax=basehisto->GetBinContent(3+i);
                bmin2=basehisto->GetBinContent((Int_t)basehisto->GetBinContent(1)+1+2+j);
                bmax2=basehisto->GetBinContent((Int_t)basehisto->GetBinContent(1)+1+3+j);
                numhisto+=1;
//              printf("%d : Min = %f Max = %f, Min2 = %f Max2 = %f\n",numhisto,bmin,bmax,bmin2,bmax2);
                if((binvalue>bmin&&binvalue<bmax) && (binvalue2>bmin2&&binvalue2<bmax2))
                        ((TH1F*)fHisto1D[SysNum]->At(numhisto))->Fill(value,wait); 
        }
      }
        
 } else {
    if(SysNum==0) {
       ostringstream error;
       error << "THistoManager::Fill1BinHisto = Error bin histograms hnum " 
             << hnum << " doesn't exist "<< endl;              
       string mess=error.str();      
       throw mess; 
    }
 }
} 




void THistoManager::ReadParamsh1f(istringstream &inParams,string type, Int_t hoffset,
                        const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo h1f

   int num=0,bin=0;
   string name,word;
   string comment="no comment";
   string space=" ";
   double min=0,max=0;

   inParams >> num >> name  >> word;
   num += hoffset;
   
   if(word=="'"){
     word=space; 
     comment=space;
     do {
       comment=comment+space+word; 
       inParams >> word;
    
     }while(word!="'"); 
   }
   if(comment.find("'")!=string::npos){
     cout << "THistoManager::ReadParameters = error in h1f " << num<< " comment format"<<endl; 
     throw string("stop! \n");
   }                 
   inParams  >>bin>>min>>max;
  

   name="h"+type+name;
   comment=comment+space+"("+type+")";
   if(SysNum==0) cout << num <<" " <<name <<" : "<<comment<<endl;
           
   if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(name.c_str(),comment.c_str(),bin,min,max); 
        ((TH1F*)fHisto1D[SysNum]->At(num))->Sumw2();
   }      


}

void THistoManager::ReadParamsh2f(istringstream &inParams,string type, Int_t hoffset,
                                const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo h2f

   int num=0,bin_x=0,bin_y=0;
   string name,word;
   string comment="no comment";
   string space=" ";
   double min_x=0,max_x=0, min_y=0,max_y=0;

   inParams >> num >> name  >> word;
   num += hoffset;
   
   if(word=="'"){
     word=space; 
     comment=space;
     do {
        comment=comment+space+word; 
        inParams >> word;
     }while(word!="'"); 
   } 
   if(comment.find("'")!=string::npos){
     cout << "THistoManager::ReadParameters = error in h2f " << num<< " comment format"<<endl; 
     throw string("stop! \n");
   }         
   inParams   >>bin_x>>min_x>>max_x>>bin_y>>min_y>>max_y; 

   name="h"+type+name;
   comment=comment+space+"("+type+")";

   if(SysNum==0) cout <<num <<" " << name <<": "<<comment <<endl;          
   if(SysNum==0||BookforSys) {
        new (fHisto2D[SysNum]->AddrAt(num)) 
                TH2F(name.c_str(),comment.c_str(),bin_x,min_x,max_x,bin_y,min_y,max_y); 
        ((TH2F*)fHisto2D[SysNum]->At(num))->Sumw2();
   }

}





void THistoManager::ReadParams1binfh1f(istringstream &inParams,string type, Int_t hoffset,
                                        const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo 1binfh1f

    cout << "Binned histo with fix bin width ..."<<endl;
    int num=0,bin=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    double min=0,max=0;
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
   inParams  >>bin>>min>>max;
    cout <<num<<" "<<name<<" "<<comment<<" "<<bin<<" "<<min<<" "<<max<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,hnum;
    Double_t binwidth,bmin,bmax,offset;
    string binvar;
    Char_t hname[100],label[200];
    inParams >> binvar >> nbbins >> binwidth >> offset; 
    cout<<" "<<binvar<<" "<< nbbins<<" "<<binwidth<<endl; 
    sprintf(hname,"h%s%sbin%s",type.c_str(),name.c_str(),binvar.c_str());
    sprintf(label,"h%s%s bin values",type.c_str(),comment.c_str());
    if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(hname,label,nbbins+2,0,nbbins+2); //histo to store min/max of bins
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(1,nbbins);
    }
    for(Int_t i=0;i<nbbins;i++) {
          bmin=offset+i*binwidth;
          bmax=offset+(i+1)*binwidth;
          hnum=num+i+1;
          sprintf(hname,"h%s%s%d",type.c_str(),name.c_str(),i);
          sprintf(label,"%s (%7.2f<%s<%7.2f) (%s)",comment.c_str(),bmin,binvar.c_str(),bmax,type.c_str());
          cout << hname <<": "<< label <<endl;          
          if(SysNum==0||BookforSys) {
                new (fHisto1D[SysNum]->AddrAt(hnum)) TH1F(hname,label,bin,min,max);
                ((TH1F*)fHisto1D[SysNum]->At(hnum))->Sumw2();      
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+2,bmin);
                if(i==nbbins-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+3,bmax);
          }
    }

}


void THistoManager::ReadParamsh1fgeneric(istringstream &inParams,string type, Int_t hoffset,
                                         list<string>& PartTags,
                                         const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo h1fgeneric

    cout << "Generic h1f histo ..."<<endl;
    int num=0,bin=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    double min=0,max=0;

    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
        comment=comment+space+word; 
        inParams >> word;
     
      }while(word!="'"); 
    }
    if(comment.find("'")!=string::npos){
      cout << "THistoManager::ReadParameters = error in h1fgeneric " << num<< " comment format"<<endl; 
      throw string("stop! \n");
    }                 
    inParams  >>bin>>min>>max;
  
//--- loop over particles to book
   Int_t nbpart=0;
   for(list<string>::iterator it=PartTags.begin();it!=PartTags.end(); it++) {
        string& parttag=(*it);
 
        string thisname="h"+type+name+parttag;
        string thiscomment=comment+parttag+space+"("+type+")";
        Int_t thisnum=num+nbpart;
        if(SysNum==0) cout << thisnum <<" " <<thisname <<" : "<<thiscomment<<endl;
            
        if(SysNum==0||BookforSys) {
                new (fHisto1D[SysNum]->AddrAt(thisnum)) TH1F(thisname.c_str(),thiscomment.c_str(),bin,min,max); 
                ((TH1F*)fHisto1D[SysNum]->At(thisnum))->Sumw2();
        } 
        nbpart+=1;
   }     


}




void THistoManager::ReadParams1binvh1f(istringstream &inParams,string type, Int_t hoffset,
                                        const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo 1binvh1f (variable width bins)

    cout << "Binned histo with variable bin width ..."<<endl;
    int num=0,bin=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    double min=0,max=0;
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
   inParams  >>bin>>min>>max;
    cout <<num<<" "<<name<<" "<<comment<<" "<<bin<<" "<<min<<" "<<max<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,hnum;
    Double_t bmin,bmax;
    string binvar;
    Char_t hname[100],label[200];
    inParams >> binvar >> nbbins;
    inParams >> bmin;
//  cout<<" "<<binvar<<" "<< nbbins<<" "<<binwidth<<endl; 
    sprintf(hname,"h%s%sbin%s",type.c_str(),name.c_str(),binvar.c_str());
    sprintf(label,"h%s%s bin values",type.c_str(),comment.c_str());
    if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(hname,label,nbbins+2,0,nbbins+2); //histo to store min/max of bins
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(1,nbbins);
    }
    for(Int_t i=0;i<nbbins;i++) {
          inParams >> bmax;
          hnum=num+i+1;
          sprintf(hname,"h%s%s%d",type.c_str(),name.c_str(),i);
          sprintf(label,"%s (%7.2f<%s<%7.2f) (%s)",comment.c_str(),bmin,binvar.c_str(),bmax,type.c_str());
          cout << hname <<": "<< label <<endl;          
          if(SysNum==0||BookforSys) {
                new (fHisto1D[SysNum]->AddrAt(hnum)) TH1F(hname,label,bin,min,max);
                ((TH1F*)fHisto1D[SysNum]->At(hnum))->Sumw2();
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+2,bmin);
                if(i==nbbins-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+3,bmax);
          }
          bmin=bmax;
    }


}


void THistoManager::ReadParams1binVarSizeh1f(istringstream &inParams,string type, Int_t hoffset,
                                             const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for bin histo 1binvbh1f (variable width bins)

    cout << "Histo with variable bin width ..."<<endl;

    // Helper vars
    int num=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    
    // Lecture du numero, du nom de l'histo et du titre
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }

    // Lecture du nombre de bins
    Int_t NInputBin;
    inParams  >>NInputBin;

    // Boucle sur le binning 
    Double_t BinningVar[NInputBin];
    Double_t iBinMin = 0;
    Double_t iBinMax = 0;
    for(Int_t iBin = 0;iBin<NInputBin;iBin++){
      inParams >> BinningVar[iBin];
    }
    iBinMin =  BinningVar[0];
    iBinMax =  BinningVar[NInputBin-1];

    for(Int_t iBin = 0;iBin<NInputBin;iBin++){
      //      printf("Bin%d = %g\n",iBin,BinningVar[iBin]);
    }

    // Make labels
    Char_t hname[100],label[200];
    sprintf(hname,"h%s%s",type.c_str(),name.c_str());
    sprintf(label,"%s",comment.c_str());
    cout <<num<<" "<<hname<<" "<<comment<<" ("<<NInputBin<<" bins from "<<iBinMin<<" to "<<iBinMax<< ")\n";

    // Store histo
    if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(hname,label,NInputBin-1,BinningVar);
        ((TH1F*)fHisto1D[SysNum]->At(num))->Sumw2();
    }

}


void THistoManager::ReadParams2binfh1f(istringstream &inParams,string type, Int_t hoffset,
                                        const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for 2bin histo 2binfh1f

    cout << "Double binned histo with fix bin width ..."<<endl;
    int num=0,bin=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    double min=0,max=0;
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
   inParams  >>bin>>min>>max;
    cout <<num<<" "<<name<<" "<<comment<<" "<<bin<<" "<<min<<" "<<max<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,nbbins2,hnum;
    Double_t binwidth,bmin,bmax,offset;
    Double_t binwidth2,bmin2,bmax2,offset2;
    string binvar,binvar2;
    Char_t hname[100],label[200];
    inParams >> binvar >> nbbins >> binwidth >> offset;
    inParams >> binvar2 >> nbbins2 >> binwidth2 >> offset2;
    cout<<" "<<binvar<<" "<< nbbins<<" "<<binwidth<<endl; 
    cout<<" "<<binvar2<<" "<< nbbins2<<" "<<binwidth2<<endl; 

    sprintf(hname,"h%s%sbin%s_%s",type.c_str(),name.c_str(),binvar.c_str(),binvar2.c_str());
    sprintf(label,"h%s%s bin values",type.c_str(),comment.c_str());
    if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(hname,label,(nbbins+nbbins2)+3,0,(nbbins+nbbins2)+3); //histo to store min/max of bins
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(1,nbbins);
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(2,nbbins2);
    }
    hnum=num;
    for(Int_t i=0;i<nbbins;i++) {
        bmin=offset+i*binwidth;
        bmax=offset+(i+1)*binwidth;
        for(Int_t j=0;j<nbbins2;j++) {
          bmin2=offset2+j*binwidth2;
          bmax2=offset2+(j+1)*binwidth2;
          hnum+=1;
          sprintf(hname,"h%s%s%d",type.c_str(),name.c_str(),hnum-num-1);
          sprintf(label,"%s (%7.2f<%s<%7.2f) (%7.2f<%s<%7.2f) (%s)"
                       ,comment.c_str(),bmin,binvar.c_str(),bmax,bmin2,binvar2.c_str(),bmax2,type.c_str());
          cout << hname <<": "<< label <<endl;          
          if(SysNum==0||BookforSys) {
                new (fHisto1D[SysNum]->AddrAt(hnum)) TH1F(hname,label,bin,min,max);
                ((TH1F*)fHisto1D[SysNum]->At(hnum))->Sumw2();
          }
        }
        if(SysNum==0||BookforSys) {
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+3,bmin);
                if(i==nbbins-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+4,bmax);
        }
    }
    for(Int_t j=0;j<nbbins2;j++) {
        bmin2=offset2+j*binwidth2;
        bmax2=offset2+(j+1)*binwidth2;
        if(SysNum==0||BookforSys) {
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(j+nbbins+1+3,bmin2);
                if(j==nbbins2-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(j+nbbins+1+4,bmax2);
        }
    }
}




void THistoManager::ReadParams2binvh1f(istringstream &inParams,string type, Int_t hoffset,
                                        const Short_t SysNum,const Bool_t BookforSys)
{
//     read params for 2bin histo 2binvh1f variable bin width

    cout << "Double binned histo with variable bin width ..."<<endl;
    int num=0,bin=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    double min=0,max=0;
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
    inParams  >>bin>>min>>max;
    cout <<num<<" "<<name<<" "<<comment<<" "<<bin<<" "<<min<<" "<<max<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,nbbins2,hnum;
    Double_t bmin,bmax;
    Double_t bmin2,bmax2;
    string binvar,binvar2;
    Char_t hname[100],label[200];
    inParams >> binvar >> nbbins;
    Float_t *bins1=new Float_t[nbbins+1];
    inParams >> bins1[0];
    for(Int_t i=0;i<nbbins;i++) inParams >> bins1[i+1];
    inParams >> binvar2 >> nbbins2 ;
    Float_t *bins2=new Float_t[nbbins2+1];
    inParams >> bins2[0];
    for(Int_t i=0;i<nbbins2;i++) inParams >> bins2[i+1];

    cout<<" "<<binvar<<" "<< nbbins<<endl; 
    cout<<" "<<binvar2<<" "<< nbbins2<<endl; 

    sprintf(hname,"h%s%sbin%s_%s",type.c_str(),name.c_str(),binvar.c_str(),binvar2.c_str());
    sprintf(label,"h%s%s bin values",type.c_str(),comment.c_str());
    if(SysNum==0||BookforSys) {
        new (fHisto1D[SysNum]->AddrAt(num)) TH1F(hname,label,(nbbins+nbbins2)+3,0,(nbbins+nbbins2)+3); //histo to store min/max of bins
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(1,nbbins);
        ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(2,nbbins2);
    }
    hnum=num;
    
    
    for(Int_t i=0;i<nbbins;i++) {
        bmin=bins1[i];
        bmax=bins1[i+1];
        for(Int_t j=0;j<nbbins2;j++) {
          bmin2=bins2[j];
          bmax2=bins2[j+1];
          hnum+=1;
          sprintf(hname,"h%s%s%d",type.c_str(),name.c_str(),hnum-num-1);
          sprintf(label,"%s (%7.2f<%s<%7.2f) (%7.2f<%s<%7.2f) (%s)"
                       ,comment.c_str(),bmin,binvar.c_str(),bmax,bmin2,binvar2.c_str(),bmax2,type.c_str());
          cout << hname <<": "<< label <<endl;          
          if(SysNum==0||BookforSys) {
                new (fHisto1D[SysNum]->AddrAt(hnum)) TH1F(hname,label,bin,min,max);
                ((TH1F*)fHisto1D[SysNum]->At(hnum))->Sumw2();
          }
        }
        if(SysNum==0||BookforSys) {
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+3,bmin);
                if(i==nbbins-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(i+4,bmax);
        }
    }
    for(Int_t j=0;j<nbbins2;j++) {
        bmin2=bins2[j];
        bmax2=bins2[j+1];
        if(SysNum==0||BookforSys) {
                ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(j+nbbins+1+3,bmin2);
                if(j==nbbins2-1) ((TH1F*)fHisto1D[SysNum]->AddrAt(num))->SetBinContent(j+nbbins+1+4,bmax2);
        }
    }

delete [] bins1;
delete [] bins2;

}







void THistoManager::ReadParams2dbin(istringstream &inParams,string type, Int_t hoffset,const Short_t SysNum)
{
//     read params for 2d bin (by B.P. used in NCCC analysis 
    cout << "2d histo with double variable bin width ..."<<endl;
  
    int num=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
    //   inParams  >>bin>>min>>max;
    name="h"+type+name;
    cout <<num<<" "<<name<<" "<<comment<<" "<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,nbbins2;
    //    Double_t bmin,bmax;
    //   Double_t bmin2,bmax2;
    string binvar,binvar2;
//   Char_t hname[100],label[200];
    inParams >> binvar >> nbbins;
    Float_t *bins1=new Float_t[nbbins+1];
    inParams >> bins1[0];
    for(Int_t i=0;i<nbbins;i++) inParams >> bins1[i+1];
    inParams >> binvar2 >> nbbins2 ;
    Float_t *bins2=new Float_t[nbbins2+1];
    inParams >> bins2[0];
    for(Int_t i=0;i<nbbins2;i++) inParams >> bins2[i+1];

    cout<<" "<<binvar<<" "<< nbbins<<endl; 
    for(Int_t i=0;i<nbbins+1;i++){
      cout << " " << bins1[i];}
    cout << endl;
    cout<<" "<<binvar2<<" "<< nbbins2<<endl; 
   for(Int_t i=0;i<nbbins2+1;i++){
      cout << " " << bins2[i];}
    cout << endl;
// now create the histogram ...
    new (fHisto2D[SysNum]->AddrAt(num)) 
             TH2F(name.c_str(),comment.c_str(),nbbins,bins1,nbbins2,bins2);   
          ((TH2F*)fHisto2D[SysNum]->At(num))->Sumw2();


delete [] bins1;
delete [] bins2;

}

void THistoManager::ReadParams3dbin(istringstream &inParams,string type, Int_t hoffset,const Short_t SysNum)
{
//     read params for 3d bin (by B.P. used in NCCC analysis 
    cout << "3d histo with double variable bin width ..."<<endl;
  
    int num=0;
    string name,word;
    string comment="no comment";
    string space=" ";
    
    inParams >> num >> name  >> word;
    num += hoffset;
    if(word=="'"){
      word=space; 
      comment=space;
      do {
         comment=comment+space+word; 
         inParams >> word;
      }while(word!="'"); 
    }
    //   inParams  >>bin>>min>>max;
    name="h"+type+name;
    cout <<num<<" "<<name<<" "<<comment<<" "<<endl;

    //--- read binned variable, number of bins, bin width
    Int_t nbbins,nbbins2,nbbins3;
    //    Double_t bmin,bmax;
    //   Double_t bmin2,bmax2;
    string binvar,binvar2,binvar3;
//   Char_t hname[100],label[200];
    inParams >> binvar >> nbbins;
    Float_t *bins1=new Float_t[nbbins+1];
    inParams >> bins1[0];
    for(Int_t i=0;i<nbbins;i++) inParams >> bins1[i+1];

    inParams >> binvar2 >> nbbins2 ;
    Float_t *bins2=new Float_t[nbbins2+1];
    inParams >> bins2[0];
    for(Int_t i=0;i<nbbins2;i++) inParams >> bins2[i+1];

    inParams >> binvar3 >> nbbins3 ;
    Float_t *bins3=new Float_t[nbbins3+1];
    inParams >> bins3[0];
    for(Int_t i=0;i<nbbins3;i++) inParams >> bins3[i+1];


    cout<<" "<<binvar<<" "<< nbbins<<endl; 
    for(Int_t i=0;i<nbbins+1;i++){
      cout << " " << bins1[i];}
    cout << endl;
    cout<<" "<<binvar2<<" "<< nbbins2<<endl; 
   for(Int_t i=0;i<nbbins2+1;i++){
      cout << " " << bins2[i];}
    cout << endl;
    cout<<" "<<binvar3<<" "<< nbbins3<<endl; 
   for(Int_t i=0;i<nbbins3+1;i++){
      cout << " " << bins3[i];}
// now create the histogram ...
    new (fHisto3D[SysNum]->AddrAt(num)) 
             TH3F(name.c_str(),comment.c_str(),nbbins,bins1,nbbins2,bins2,nbbins3,bins3);   
          ((TH3F*)fHisto3D[SysNum]->At(num))->Sumw2();


delete [] bins1;
delete [] bins2;
delete [] bins3;

}


void THistoManager::Scale1dHisto(Int_t hnum,Float_t ScaleParam,const Short_t SysNum)
{
// Scale 1D histogram

  // Get histogram number hnum
  TH1F *Histo = (TH1F *)fHisto1D[SysNum]->At(hnum);
  Histo->Scale(ScaleParam);
}

---