[u/mrichter/AliRoot.git] / PWG4 / TwoPartCorr / mergeCorr.C

// $Id:$

#if !defined(__CINT__) || defined(__MAKECINT__)
#include <vector>
#include <map>
#include <Riostream.h>
#include <TChain.h>
#include <TClonesArray.h>
#include <TError.h>
#include <TFile.h>
#include <TDirectory.h>
#include <TMath.h>
#include <TRandom3.h>
#include <TSystem.h>
#include "MyClasses.cxx"
#endif

//#define DUPLICATECHECK
#ifdef DUPLICATECHECK
//#define DUPLICATEVERBOSE
#endif

struct EvInfo
{
  Int_t         fRun;
  Short_t       fArrId;
  UInt_t        fEntry; 
};


void mergeCorr(Int_t run_from=-1, 
               Int_t run_to=-1,
               const char *outPrefix = "merged",
               const char *inFileNames = "res/output_*.root");

void mergeCorr(Int_t nEvents,
               const char *outFileName,
               const char *inFileNames);

//-----------------------------------------------------------------------------------------------------

void mergeCorr(Int_t run_from,
               Int_t run_to,
               const char *outPrefix,
               const char *inFileNames)
{
  TChain *c = new TChain("MyTree");
  c->Add(inFileNames);

  std::map<Int_t,Int_t> goodruns;

  TObjArray *filenames = c->GetListOfFiles();
  Int_t nfiles = filenames->GetEntries();
  for(Int_t i=0; i<nfiles; ++i) {
    TString fname(filenames->At(i)->GetTitle());
    TFile *file = TFile::Open(fname);
    if (!file)
      continue;
    TList *output = dynamic_cast<TList*>(file->Get("output"));
    TTree *tree = dynamic_cast<TTree*>(output->FindObject("MyTree"));
    Int_t nent = tree->GetEntries();
    if (nent<1) {
      delete tree;
      file->Close();
      delete file;
      continue;
    }

    MyHeader *header = 0;
    TBranch *br = tree->GetBranch("header");
    br->SetAddress(&header);

    std::vector<Int_t> multruns;

    Int_t runno = -1;
    for (Int_t ev=0;ev<nent;++ev) {
      br->GetEntry(ev);
      if (header->fRun==0) {
        continue;
      }
      if (runno<0) 
        runno = header->fRun;
      else if (runno!=header->fRun) {
        Bool_t found = 0;
        for (UInt_t j=0;j<multruns.size();++j) {
          if (multruns.at(j)==header->fRun) {
            found = 1;
            break;
          }
        }
        if (!found)
          multruns.push_back(header->fRun);
      }
    }
    if (multruns.size()<=0) {
      if ((runno>run_from) && (run_to<0 || runno<run_to)) {
        cout << "Found run " << runno << " associated to " << fname << endl;
        goodruns.insert(pair<Int_t,Int_t>(runno,0));
        TString dir(gSystem->DirName(inFileNames));
        dir+="/runs";
        TString link(Form("%s/%s",gSystem->pwd(),fname.Data()));
        TString rundir(Form("%s/%d", dir.Data(), runno));
        TString cmd(Form("mkdir -p %s && cd %s && ln -sf %s", 
                         rundir.Data(), rundir.Data(), link.Data()));
        gSystem->Exec(cmd);
      }
    } else {
      cout << "Multiple runs (" << runno;
      for (UInt_t j=0;j<multruns.size();++j)
        cout << ", " << multruns.at(j);
      cout << ") found in " << fname << endl;
      cout << " Cannot deal with this case yet!!!" << endl;
    }
    delete tree;
    file->Close();
    delete file;
  }

  for (map<Int_t, Int_t>::iterator it = goodruns.begin();it!=goodruns.end();++it) {
    TString infiles(Form("%s/runs/%d/*.root",gSystem->DirName(inFileNames),it->first));
    TString outdir(Form("%s/mergedruns",gSystem->DirName(inFileNames)));
    TString outFileName(Form("%s/%s_run%d.root",outdir.Data(),outPrefix,it->first));
    gSystem->Exec(Form("mkdir -p %s",outdir.Data()));
    gSystem->Exec(Form("root -b -q mergeCorr.C+'(-1,\"%s\",\"%s\")'",
                       outFileName.Data(),infiles.Data()));
  }
}

//-----------------------------------------------------------------------------------------------------

void mergeCorr(Int_t nEvents,
               const char *outFileName,
               const char *inFileNames)
{
  TChain *c = new TChain("MyTree");
  c->Add(inFileNames);

  map<ULong64_t, EvInfo*> einfos;
#ifdef DUPLICATEVERBOSE
  map<ULong64_t, MyHeader*> eheaders;
#endif

  Int_t totnent = 0;
  Int_t totnsus = 0;
#ifdef DUPLICATECHECK
  Int_t totncan = 0;
#ifdef DUPLICATEVERBOSE
  Int_t totndup = 0;
#endif
#endif

  TObjArray *filenames = c->GetListOfFiles();
  Int_t nfiles = filenames->GetEntries();
  for(Int_t i=0; i<nfiles; ++i) {
    TString fname(filenames->At(i)->GetTitle());
    TFile *file = TFile::Open(fname);
    if (!file)
      continue;
    TList *output = dynamic_cast<TList*>(file->Get("output"));
    TTree *tree = dynamic_cast<TTree*>(output->FindObject("MyTree"));
    Int_t nent = tree->GetEntries();
    if (nent<1) {
      delete tree;
      file->Close();
      delete file;
      continue;
    }

    totnent += nent;
    cout << "Found "<< nent << " entries in " << fname << endl;

    MyHeader *header = 0;
    TBranch *br = tree->GetBranch("header");
    br->SetAddress(&header);

    for (Int_t ev=0;ev<nent;++ev) {
      br->GetEntry(ev);
      if (header->fPeriod!=0) // fix take period from runinfo block
        cout << "--> Non-zero period found " << header->fRun << " " 
             << header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
      if (header->fRun==0) {
#ifdef DUPLICATEVERBOSE
        cout << "--> Suspicious found for " << header->fRun << " " 
             << header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
#endif
        ++totnsus;
        continue;
      }

      ULong64_t evtid = header->fTime;//header->GetEventId();
#ifdef DUPLICATECHECK
      if (einfos.find(evtid)!=einfos.end()) {
#ifdef DUPLICATEVERBOSE        
        cout << "--> Potential duplicate found for " << header->fRun << " " 
             << header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
#endif
      ++totncan;
#ifdef DUPLICATEVERBOSE        
        map<ULong64_t, MyHeader*>::iterator it = eheaders.find(evtid);
        if ( (it->second->fRun==header->fRun) && 
             (it->second->fNTracks==header->fNTracks) &&
             (it->second->fNSelTracks==header->fNSelTracks) &&
             (it->second->fNTracklets==header->fNTracklets) &&
             (it->second->fVz==header->fVz)) {
          cout << "--> Duplicate found:" << endl;
          cout << " Run: " << it->second->fRun << " vs " << header->fRun << endl; 
          cout << " Ntracks: " << it->second->fNTracks << " vs " << header->fNTracks << endl; 
          cout << " Nseltracks: " << it->second->fNSelTracks << " vs " << header->fNSelTracks << endl; 
          cout << " Ntracklets: " << it->second->fNTracklets << " vs " << header->fNTracklets << endl; 
          cout << " Vz: " << it->second->fVz << " vs " << header->fVz << endl;
          ++totndup;
        }
#endif
        continue;
      }
#ifdef DUPLICATEVERBOSE        
      MyHeader *cheader = new MyHeader(*header);
      eheaders.insert( pair<ULong64_t, MyHeader*>(evtid,cheader) );
#endif
#endif
      EvInfo *cinfo = new EvInfo;
      cinfo->fRun = header->fRun;
      cinfo->fArrId = i;
      cinfo->fEntry = ev;
      einfos.insert( pair<ULong64_t, EvInfo*>(evtid,cinfo) );
    }
    delete tree;
    file->Close();
    delete file;
    if ((nEvents>0) && (einfos.size()>=(UInt_t)nEvents))
      break;
  }
  cout << "Found total: " << totnent << endl;
  cout << "Found suspicious: " << totnsus << endl;
#ifdef DUPLICATECHECK
  cout << "Found potential duplicates: " << totncan << endl;
#ifdef DUPLICATEVERBOSE        
  cout << "Found duplicates: " << totndup << endl;
#endif
#endif

  // sort output tree
  std::vector<TFile*> files(nfiles);
  std::vector<TTree*> trees(nfiles);
  std::vector<MyHeader*> headers(nfiles);
  std::vector<TClonesArray*> partas(nfiles);

  MyHeader *newHeader = new MyHeader;
  if (TClass::GetClass("MyPart"))
    TClass::GetClass("MyPart")->IgnoreTObjectStreamer();
  TClonesArray *newParts = new TClonesArray("MyPart",10000);

  TFile *newFile = TFile::Open(outFileName,"recreate","Merged/Sorted MyTree created by mergeCorr.C",5);
  TDirectory::TContext cd(newFile);
  TTree *newTree = new TTree("MyTree", "MyTree created by MyAliTask.cxx");
  newTree->Branch("header",&newHeader, 32*1024, 99);
  newTree->Branch("parts",&newParts, 256*1024, 99);
  Int_t newEntries = 0;
  for (map<ULong64_t, EvInfo*>::iterator it = einfos.begin();it!=einfos.end();++it) {
    Short_t id = it->second->fArrId;
    UInt_t entry = it->second->fEntry;
    if (!trees.at(id)) { // connect new tree
      TString fname(filenames->At(id)->GetTitle());
      TFile *file = TFile::Open(fname);
      TDirectory::TContext cd2(newFile,file);
      files[id] = file;
      TList *output = dynamic_cast<TList*>(file->Get("output"));
      TTree *tree = dynamic_cast<TTree*>(output->FindObject("MyTree"));
      trees[id] = tree;
      tree->SetBranchAddress("header",&headers[id]);
      tree->SetBranchAddress("parts",&partas[id]);
    }
    files[id]->cd();
    trees.at(id)->GetEntry(entry);
    *newHeader = *(headers[id]);
    newParts->Clear();
    Int_t nparts = partas[id]->GetEntries();
    for (Int_t p = 0; p < nparts; ++p) {
      MyPart *orig = dynamic_cast<MyPart*>(partas[id]->At(p));
      new((*newParts)[p]) MyPart(*orig);
    }
    newFile->cd();
    newTree->Fill();
    ++newEntries;
    delete it->second;
    if ((nEvents>0) && (newEntries>=nEvents))
      break;
  }

  newTree->Write();
  newFile->Close();
  delete newFile;

  for(Int_t i=0;i<nfiles;++i) {
    delete trees.at(i);
    delete files.at(i);
  }
}
Commit	Line	Data
828cd110	1	// $Id:$
	2
	3	#if !defined(__CINT__) \|\| defined(__MAKECINT__)
	4	#include <vector>
	5	#include <map>
	6	#include <Riostream.h>
	7	#include <TChain.h>
	8	#include <TClonesArray.h>
	9	#include <TError.h>
	10	#include <TFile.h>
	11	#include <TDirectory.h>
	12	#include <TMath.h>
	13	#include <TRandom3.h>
	14	#include <TSystem.h>
	15	#include "MyClasses.cxx"
	16	#endif
	17
	18	//#define DUPLICATECHECK
	19	#ifdef DUPLICATECHECK
	20	//#define DUPLICATEVERBOSE
	21	#endif
	22
	23	struct EvInfo
	24	{
	25	Int_t fRun;
	26	Short_t fArrId;
	27	UInt_t fEntry;
	28	};
	29
	30
	31	void mergeCorr(Int_t run_from=-1,
	32	Int_t run_to=-1,
	33	const char *outPrefix = "merged",
	34	const char inFileNames = "res/output_.root");
	35
	36	void mergeCorr(Int_t nEvents,
	37	const char *outFileName,
	38	const char *inFileNames);
	39
	40	//-----------------------------------------------------------------------------------------------------
	41
	42	void mergeCorr(Int_t run_from,
	43	Int_t run_to,
	44	const char *outPrefix,
	45	const char *inFileNames)
	46	{
	47	TChain *c = new TChain("MyTree");
	48	c->Add(inFileNames);
	49
	50	std::map<Int_t,Int_t> goodruns;
	51
	52	TObjArray *filenames = c->GetListOfFiles();
	53	Int_t nfiles = filenames->GetEntries();
	54	for(Int_t i=0; i<nfiles; ++i) {
	55	TString fname(filenames->At(i)->GetTitle());
	56	TFile *file = TFile::Open(fname);
	57	if (!file)
	58	continue;
	59	TList output = dynamic_cast<TList>(file->Get("output"));
	60	TTree tree = dynamic_cast<TTree>(output->FindObject("MyTree"));
	61	Int_t nent = tree->GetEntries();
	62	if (nent<1) {
497ecc12	63	delete tree;
828cd110	64	file->Close();
	65	delete file;
	66	continue;
	67	}
	68
	69	MyHeader *header = 0;
	70	TBranch *br = tree->GetBranch("header");
	71	br->SetAddress(&header);
	72
	73	std::vector<Int_t> multruns;
	74
	75	Int_t runno = -1;
	76	for (Int_t ev=0;ev<nent;++ev) {
	77	br->GetEntry(ev);
	78	if (header->fRun==0) {
	79	continue;
	80	}
	81	if (runno<0)
	82	runno = header->fRun;
	83	else if (runno!=header->fRun) {
	84	Bool_t found = 0;
	85	for (UInt_t j=0;j<multruns.size();++j) {
	86	if (multruns.at(j)==header->fRun) {
	87	found = 1;
	88	break;
	89	}
	90	}
	91	if (!found)
	92	multruns.push_back(header->fRun);
	93	}
	94	}
	95	if (multruns.size()<=0) {
	96	if ((runno>run_from) && (run_to<0 \|\| runno<run_to)) {
	97	cout << "Found run " << runno << " associated to " << fname << endl;
	98	goodruns.insert(pair<Int_t,Int_t>(runno,0));
	99	TString dir(gSystem->DirName(inFileNames));
	100	dir+="/runs";
	101	TString link(Form("%s/%s",gSystem->pwd(),fname.Data()));
	102	TString rundir(Form("%s/%d", dir.Data(), runno));
	103	TString cmd(Form("mkdir -p %s && cd %s && ln -sf %s",
	104	rundir.Data(), rundir.Data(), link.Data()));
	105	gSystem->Exec(cmd);
	106	}
	107	} else {
	108	cout << "Multiple runs (" << runno;
	109	for (UInt_t j=0;j<multruns.size();++j)
	110	cout << ", " << multruns.at(j);
	111	cout << ") found in " << fname << endl;
	112	cout << " Cannot deal with this case yet!!!" << endl;
	113	}
497ecc12	114	delete tree;
828cd110	115	file->Close();
	116	delete file;
	117	}
	118
	119	for (map<Int_t, Int_t>::iterator it = goodruns.begin();it!=goodruns.end();++it) {
	120	TString infiles(Form("%s/runs/%d/*.root",gSystem->DirName(inFileNames),it->first));
	121	TString outdir(Form("%s/mergedruns",gSystem->DirName(inFileNames)));
	122	TString outFileName(Form("%s/%s_run%d.root",outdir.Data(),outPrefix,it->first));
	123	gSystem->Exec(Form("mkdir -p %s",outdir.Data()));
	124	gSystem->Exec(Form("root -b -q mergeCorr.C+'(-1,\"%s\",\"%s\")'",
	125	outFileName.Data(),infiles.Data()));
	126	}
	127	}
	128
	129	//-----------------------------------------------------------------------------------------------------
	130
	131	void mergeCorr(Int_t nEvents,
	132	const char *outFileName,
	133	const char *inFileNames)
	134	{
	135	TChain *c = new TChain("MyTree");
	136	c->Add(inFileNames);
	137
	138	map<ULong64_t, EvInfo*> einfos;
	139	#ifdef DUPLICATEVERBOSE
	140	map<ULong64_t, MyHeader*> eheaders;
	141	#endif
	142
	143	Int_t totnent = 0;
	144	Int_t totnsus = 0;
	145	#ifdef DUPLICATECHECK
	146	Int_t totncan = 0;
	147	#ifdef DUPLICATEVERBOSE
	148	Int_t totndup = 0;
	149	#endif
	150	#endif
	151
	152	TObjArray *filenames = c->GetListOfFiles();
	153	Int_t nfiles = filenames->GetEntries();
	154	for(Int_t i=0; i<nfiles; ++i) {
	155	TString fname(filenames->At(i)->GetTitle());
	156	TFile *file = TFile::Open(fname);
	157	if (!file)
	158	continue;
	159	TList output = dynamic_cast<TList>(file->Get("output"));
	160	TTree tree = dynamic_cast<TTree>(output->FindObject("MyTree"));
	161	Int_t nent = tree->GetEntries();
	162	if (nent<1) {
497ecc12	163	delete tree;
828cd110	164	file->Close();
	165	delete file;
	166	continue;
	167	}
	168
	169	totnent += nent;
	170	cout << "Found "<< nent << " entries in " << fname << endl;
	171
	172	MyHeader *header = 0;
	173	TBranch *br = tree->GetBranch("header");
	174	br->SetAddress(&header);
	175
	176	for (Int_t ev=0;ev<nent;++ev) {
	177	br->GetEntry(ev);
	178	if (header->fPeriod!=0) // fix take period from runinfo block
	179	cout << "--> Non-zero period found " << header->fRun << " "
	180	<< header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
	181	if (header->fRun==0) {
	182	#ifdef DUPLICATEVERBOSE
	183	cout << "--> Suspicious found for " << header->fRun << " "
	184	<< header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
	185	#endif
	186	++totnsus;
	187	continue;
	188	}
	189
	190	ULong64_t evtid = header->fTime;//header->GetEventId();
	191	#ifdef DUPLICATECHECK
	192	if (einfos.find(evtid)!=einfos.end()) {
	193	#ifdef DUPLICATEVERBOSE
	194	cout << "--> Potential duplicate found for " << header->fRun << " "
	195	<< header->fTime << " " << header->fPeriod << " " << header->fBx << endl;
	196	#endif
	197	++totncan;
	198	#ifdef DUPLICATEVERBOSE
	199	map<ULong64_t, MyHeader*>::iterator it = eheaders.find(evtid);
	200	if ( (it->second->fRun==header->fRun) &&
	201	(it->second->fNTracks==header->fNTracks) &&
	202	(it->second->fNSelTracks==header->fNSelTracks) &&
	203	(it->second->fNTracklets==header->fNTracklets) &&
	204	(it->second->fVz==header->fVz)) {
	205	cout << "--> Duplicate found:" << endl;
	206	cout << " Run: " << it->second->fRun << " vs " << header->fRun << endl;
	207	cout << " Ntracks: " << it->second->fNTracks << " vs " << header->fNTracks << endl;
	208	cout << " Nseltracks: " << it->second->fNSelTracks << " vs " << header->fNSelTracks << endl;
	209	cout << " Ntracklets: " << it->second->fNTracklets << " vs " << header->fNTracklets << endl;
	210	cout << " Vz: " << it->second->fVz << " vs " << header->fVz << endl;
	211	++totndup;
	212	}
	213	#endif
	214	continue;
	215	}
	216	#ifdef DUPLICATEVERBOSE
	217	MyHeader cheader = new MyHeader(header);
	218	eheaders.insert( pair<ULong64_t, MyHeader*>(evtid,cheader) );
	219	#endif
	220	#endif
	221	EvInfo *cinfo = new EvInfo;
	222	cinfo->fRun = header->fRun;
	223	cinfo->fArrId = i;
	224	cinfo->fEntry = ev;
	225	einfos.insert( pair<ULong64_t, EvInfo*>(evtid,cinfo) );
	226	}
497ecc12	227	delete tree;
828cd110	228	file->Close();
	229	delete file;
	230	if ((nEvents>0) && (einfos.size()>=(UInt_t)nEvents))
	231	break;
	232	}
	233	cout << "Found total: " << totnent << endl;
	234	cout << "Found suspicious: " << totnsus << endl;
	235	#ifdef DUPLICATECHECK
	236	cout << "Found potential duplicates: " << totncan << endl;
	237	#ifdef DUPLICATEVERBOSE
	238	cout << "Found duplicates: " << totndup << endl;
	239	#endif
	240	#endif
	241
	242	// sort output tree
	243	std::vector<TFile*> files(nfiles);
	244	std::vector<TTree*> trees(nfiles);
	245	std::vector<MyHeader*> headers(nfiles);
	246	std::vector<TClonesArray*> partas(nfiles);
	247
	248	MyHeader *newHeader = new MyHeader;
	249	if (TClass::GetClass("MyPart"))
	250	TClass::GetClass("MyPart")->IgnoreTObjectStreamer();
	251	TClonesArray *newParts = new TClonesArray("MyPart",10000);
	252
	253	TFile *newFile = TFile::Open(outFileName,"recreate","Merged/Sorted MyTree created by mergeCorr.C",5);
	254	TDirectory::TContext cd(newFile);
	255	TTree *newTree = new TTree("MyTree", "MyTree created by MyAliTask.cxx");
	256	newTree->Branch("header",&newHeader, 32*1024, 99);
	257	newTree->Branch("parts",&newParts, 256*1024, 99);
	258	Int_t newEntries = 0;
	259	for (map<ULong64_t, EvInfo*>::iterator it = einfos.begin();it!=einfos.end();++it) {
	260	Short_t id = it->second->fArrId;
	261	UInt_t entry = it->second->fEntry;
	262	if (!trees.at(id)) { // connect new tree
	263	TString fname(filenames->At(id)->GetTitle());
	264	TFile *file = TFile::Open(fname);
	265	TDirectory::TContext cd2(newFile,file);
	266	files[id] = file;
	267	TList output = dynamic_cast<TList>(file->Get("output"));
	268	TTree tree = dynamic_cast<TTree>(output->FindObject("MyTree"));
	269	trees[id] = tree;
	270	tree->SetBranchAddress("header",&headers[id]);
	271	tree->SetBranchAddress("parts",&partas[id]);
	272	}
	273	files[id]->cd();
	274	trees.at(id)->GetEntry(entry);
	275	newHeader = (headers[id]);
	276	newParts->Clear();
	277	Int_t nparts = partas[id]->GetEntries();
	278	for (Int_t p = 0; p < nparts; ++p) {
	279	MyPart orig = dynamic_cast<MyPart>(partas[id]->At(p));
	280	new((newParts)[p]) MyPart(orig);
	281	}
	282	newFile->cd();
	283	newTree->Fill();
	284	++newEntries;
	285	delete it->second;
	286	if ((nEvents>0) && (newEntries>=nEvents))
	287	break;
	288	}
	289
	290	newTree->Write();
	291	newFile->Close();
292	delete newFile;
293
497ecc12	294	for(Int_t i=0;i<nfiles;++i) {
497ecc12	295	delete trees.at(i);
828cd110	296	delete files.at(i);
497ecc12	297	}
828cd110	298	}