MadChain.cxx

Go to the documentation of this file.

#ifndef madchain_cxx
#define madchain_cxx
#include <iostream>
//root includes
#include "Mad/MadChain.h"

//ClassImp(MadChain);

MadChain::MadChain(TChain *chainSR,TChain *chainMC,
                   TChain *chainTH,TChain *chainEM, bool build_lookup)
  : fLookupBuilt(false)
{
  
  Init(chainSR,chainMC,chainTH,chainEM);
  if(build_lookup) BuildLookup();
}

MadChain::~MadChain()
{
  if(!fChainSR) return;
  delete fChainSR->GetCurrentFile();
  if(fChainEM) delete fChainEM->GetCurrentFile();
  if(!fChainMC) return;
  delete fChainMC->GetCurrentFile();
  if(!fChainTH) return;
  delete fChainTH->GetCurrentFile();  
}

Int_t MadChain::GetEntryNumber(Int_t run,Int_t snarl)
{
  if(!fLookupBuilt) return fChainSR->GetEntryNumberWithIndex(run,snarl);
  else{
    std::pair<Int_t,Int_t> p(run,snarl);
    std::map< std::pair<Int_t,Int_t> , Long64_t>::const_iterator it 
      = fLookUp.find(p);
    if(it==fLookUp.end()) return -1;
    else return it->second;
  }
}

Int_t MadChain::Get(Int_t entry)
{
  
  //estimate MC,TH entry numbers
  if(entry<fCurrentSR||fCurrentSR==-1) {
    fCurrentMC=entry;
    fCurrentTH=entry;
    fCurrentEM=entry;
  }
  else if(entry>fCurrentSR) {
    fCurrentMC+=(entry-fCurrentSR);
    fCurrentTH+=(entry-fCurrentSR);
    fCurrentEM+=(entry-fCurrentSR);
  }
  
  // Read contents of entry.
  if (!fChainSR) return 0;
  Int_t status = fChainSR->GetEntry(entry);
  fCurrentSR = entry;

  if(isEM) {
    fChainEM->GetEntry(fCurrentEM);
    int xtra = 0;  //just in case tree entries are not correlated
    while(Record->GetHeader().GetSnarl()!=emRecord->GetHeader().GetSnarl()){
    xtra+=1;
    fChainEM->GetEntry(fCurrentEM+xtra);
    }
    fCurrentEM+=xtra;
  }
  
  if(isMC) {
    fChainMC->GetEntry(fCurrentMC);
    int xtra = 0;  //just in case tree entries are not correlated
    while(Record->GetHeader().GetSnarl()!=mcRecord->GetHeader().GetSnarl()){
    xtra+=1;
    fChainMC->GetEntry(fCurrentMC+xtra);     
    }
    fCurrentMC+=xtra;
  }
  
  if(isTH) {
    fChainTH->GetEntry(fCurrentTH);
    int xtra = 0;  //just in case tree entries are not correlated
    while(Record->GetHeader().GetSnarl()!=thRecord->GetHeader().GetSnarl()){
    xtra+=1;
    fChainTH->GetEntry(fCurrentTH+xtra);      
    }
    fCurrentTH+=xtra;
  }
  return status;
}

void MadChain::Init(TChain *chainSR,TChain *chainMC,
                    TChain *chainTH,TChain *chainEM)
{

  Zero();

  if(chainSR==0) {
    std::cerr << "No file added" << std::endl;
    return;
  }
  fChainSR = chainSR;
  if(strcmp(fChainSR->GetName(),"NtpSR")==0){
    fChainSR->SetBranchAddress("NtpSRRecord",&Record);
  }
  else if(strcmp(fChainSR->GetName(),"NtpSt")==0){
    isST = true;
    fChainSR->SetBranchAddress("NtpStRecord",&stRecord);
  }
  fCurrentSR = -1;
  Nentries = Int_t(fChainSR->GetEntries());
  //chainSR->BuildIndex("fHeader.fRun","fHeader.fSnarl");

  if(chainEM==0) {
    isEM = false;    
  }
  else {
    fChainEM = chainEM;
    fChainEM->SetBranchAddress("NtpEMRecord",&emRecord);
    fCurrentEM = -1;    
  }

  if(chainMC==0) {
    isMC = false;
    isTH = false;
    return;
  }
  fChainMC = chainMC;
  fChainMC->SetBranchAddress("NtpMCRecord",&mcRecord);
  fCurrentMC = -1;

  if(chainTH==0) {
    isTH = false;
    return;
  }
  fChainTH = chainTH;
  fChainTH->SetBranchAddress("NtpTHRecord",&thRecord);
  fCurrentTH = -1;

}

void MadChain::Show(Int_t entry)
{
// Print contents of entry.
// If entry is not specified, print current entry
   if (!fChainSR) return;
   fChainSR->Show(entry);
   if (fChainEM) fChainEM->Show(entry);
   if (!fChainMC) return;
   fChainMC->Show(entry);
   if (!fChainTH) return;
   fChainTH->Show(entry);
}

void MadChain::Zero(){
    
  Record = 0;
  stRecord = 0;
  emRecord = 0;
  mcRecord = 0;
  thRecord = 0;
  
  isMC = true;
  isTH = true;
  isEM = true;
  isST = false;
  Nentries = -1;
  
  fChainSR = 0;
  fChainMC = 0;
  fChainTH = 0;
  fChainEM = 0;
  
  fCurrentSR = 0;
  fCurrentMC = 0;
  fCurrentTH = 0;
  fCurrentEM = 0;

}

void MadChain::BuildLookup(){
  // fill fLookUp with run,snarl of each entry in the chain
  // do it as efficiently as possible by only reading the header

  if(!fChainSR) return;

  /*
  RecRecordImp<RecCandHeader>* recrecord = 0;
  //  fChainSR->SetBranchAddress(fChainSR->GetListOfBranches()->At(0)->GetName(),
  //                         &recrecord);
  fChainSR->SetBranchAddress("NtpStRecord",
                             &recrecord);
  fChainSR->SetBranchStatus("*",0);
  fChainSR->SetBranchStatus("fHeader.*",1);
  fChainSR->SetBranchAddress("NtpStRecord",
                             &recrecord);
  */
  /*
  fChainSR->SetBranchStatus("fHeader.fRun",1);  
  fChainSR->SetBranchStatus("fHeader.fSnarl",1);  
  Int_t run, snarl;
  if(isST) {
    fChainSR->SetBranchAddress("NtpStRecord",&stRecord);
    fChainSR->SetBranchAddress("fHeader.fRun",&run);
    fChainSR->SetBranchAddress("fHeader.fSnarl",&snarl);
  }
  else{
    fChainSR->SetBranchAddress("NtpSRRecord",&Record);
    fChainSR->SetBranchAddress("fHeader.fRun",&run);
    fChainSR->SetBranchAddress("fHeader.fSnarl",&snarl);
  }
  */


  Long64_t N = fChainSR->GetEntries();
  Int_t nbytes=0;
  std::cout<<"MadChain building lookup table"<<std::endl;
  //  Int_t current_tree=-1; //fChainSR->GetTreeNumber();
  for(Long64_t i = 0; i<N; i++){
    std::pair<Int_t,Int_t> p;
    /* // this trash doesn't work
    // mock up what MakeClass generated code does
    Long64_t centry = fChainSR->LoadTree(i);
    if(fChainSR->GetTreeNumber() != current_tree){
      current_tree=fChainSR->GetTreeNumber();
      fChainSR->SetBranchStatus("*",0);
      //      fChainSR->SetBranchStatus("fHeader.*",1);
      //      fChainSR->SetBranchAddress("NtpStRecord", &recrecord);
    }    
    nbytes+=fChainSR->GetEntry(centry);
    */
    
    nbytes+=fChainSR->GetEntry(i);
    
    if(isST) {
      p.first = stRecord->GetHeader().GetRun();
      p.second = stRecord->GetHeader().GetSnarl();
    }
    else{
      p.first = Record->GetHeader().GetRun();
      p.second = Record->GetHeader().GetSnarl();
    }

    //    p.first=run; p.second=snarl;
    //    if(i%1000==0) recrecord->GetHeader().Print();
    // p.first = recrecord->GetHeader().GetRun();
    // p.second = recrecord->GetHeader().GetSnarl();
 
    fLookUp.insert( std::pair< std::pair<Int_t,Int_t>, Long64_t>(p, i) );
    if(i%10000==0){
      std::cout<<"progress: "<<i<<" events and ";
      std::cout.precision(3);
      std::cout<<(nbytes/1024.0)<<" kb read \r";
      std::cout.flush();
    }
    //    recrecord->Clear();
  }
  std::cout<<"Done                                                "<<std::endl;
  std::cout<<"lookup table has "<<fLookUp.size()<<" entries"<<std::endl;

  /*
  std::map< std::pair<Int_t,Int_t> , Long64_t>::const_iterator it
    =fLookUp.begin();
  int cntr=0;
  for(; it!=fLookUp.end(); it++){
    if(cntr%1000==0){      
      std::cout<<cntr<<" "<<it->first.first<<" "
               <<it->first.second<<" "<<it->second<<std::endl;
    }
    cntr++;
  }
  */
  fChainSR->SetBranchStatus("*",1);
  if(isST) {
    fChainSR->SetBranchAddress("NtpStRecord",&stRecord);
  }
  else{
    fChainSR->SetBranchAddress("NtpSRRecord",&Record);
  }


  fLookupBuilt=true;
}

#endif // #ifdef madchain_cxx

---