aux/mask_profiler/mask_profiler.cpp

#include <iostream>
#include <fstream>

#include <TBenchmark.h>
#include <TMap.h>
#include <TList.h>
#include <TString.h>
#include <TObjString.h>
#include <TH1D.h>
#include <TH2D.h>
#include <TSystemDirectory.h>
#include <TSystemFile.h>
#include <TRegexp.h>

#include <PamCut/PamCutEnv.h>
#include <VKmask.h>
#include <TimeVec.h>


//Declaration of all functions
Int_t HandleInputPar(Int_t argc, char **argv, TMap* optmap); 
Bool_t IsOption(TList* optlist, TString word);
// command line parser
void PrintOptions(TMap* optmap);
//Making TH2D with LT profiles
TH2D* MakeLTmaskHisto(TObject* path_to_rig, Int_t nmasks, TimeVec* tm_vec, Int_t bun);
//Make TimeVector from textfile
TimeVec * MakeTimeVec( TObject* path_to_textfile);
Int_t FindMaskIDinList(TList *lmsk, VKmask* msk);


//general loop. Returns List of lists to store in root-files
TList * DoMasksProfile( TObject* inp_path_o, TObject* inp_pattern, TObject* path_to_rig, TimeVec* tm_vec, Double_t cut_f, Double_t b_min, Double_t b_max );

void SaveResult( TObject* out_path_o, TList* out_list);

Int_t main(Int_t argc, char **argv){
  // Arguments are following:
  //
  // --inp_dir /// input directory with root files
  // --rigbin_file /// file with left edges of rigidity bins
  // --cutof_factor /// factor to cut in SVL. Typycal value is 1.3
  // --babs_min /// minimum value of magnetic field. 0.26 to exclude SAA
  // --babs_max /// maximum value of m.f. usually 1000000 to include everything
  // --bunch_file /// textfile with bunches with start-stop values

  //Default values
  Double_t cut_f = 1.3;
  Double_t babs_min = 0.26;
  Double_t babs_max = 10000000.;
  

  TMap* optmap = new TMap(); //map with options
  HandleInputPar(argc,argv, optmap);
  PrintOptions(optmap);


  TBenchmark* bench = new TBenchmark();
  bench->Start("client");


  //Making TimeVec

  TimeVec * tm_vec = MakeTimeVec( optmap->GetValue("--bunch_file") ); 
  tm_vec->Print();


  TList *lout =  DoMasksProfile( optmap->GetValue("--inp_dir"), 
                                 optmap->GetValue("--inp_pattern"), 
                                 optmap->GetValue("--rigbin_file"),
                                 tm_vec, cut_f, babs_min, babs_max );


  SaveResult(  optmap->GetValue("--out_dir"), lout);

  // report a benchmark message to the server
  bench->Stop("client");
  TString benchresult = " job";
  benchresult.Remove(0, 1);
  benchresult += " finished (";
  TDatime d;
  benchresult += d.AsSQLString();
  benchresult += " | cpu time = ";
  benchresult += floor((Double_t) bench->GetCpuTime("client"));
  benchresult += " s = ";
  benchresult += floor(((Double_t)bench->GetCpuTime("client"))/((Double_t) bench->GetRealTime("client"))*100.);
  benchresult += " % real time)";
  cout<<benchresult.Data()<<endl;
  delete bench;


  return 0;
}


// Command line parser

Int_t HandleInputPar(Int_t argc, char **argv, TMap* optmap){


  TList* allowed_opt = new TList();
  TString  optname = "--inp_dir"; // path with pattern to input files to make chain
  allowed_opt->Add(new TObjString(optname));
  optname = "--inp_pattern"; // pattern of files to make chain
  allowed_opt->Add(new TObjString(optname));
  optname="--out_dir"; //output path to store root-files 
  allowed_opt->Add(new TObjString(optname));
  optname = "--rigbin_file"; // full path to textfile with rigidity bins
  allowed_opt->Add(new TObjString(optname));
  optname = "--cutoff_factor"; // factor for SVL. Typ. is 1.3
  allowed_opt->Add(new TObjString(optname));
  optname = "--babs_min"; // cut to remove SAA and polar regions. Typ 0.26
  allowed_opt->Add(new TObjString(optname));
  optname = "--babs_max"; // 
  allowed_opt->Add(new TObjString(optname));
  optname = "--bunch_file"; // file with start-stop abstime
  allowed_opt->Add(new TObjString(optname));


  TString arg=""; 


  if(argc>1){
    if(!strcmp(argv[1], "-h") || !strcmp(argv[1], "--help") ){
      return(1);
    };
    
    for (Int_t i = 1; i < argc; i++){
      if(IsOption(allowed_opt,argv[i])){
        TObjString* option =(TObjString*)allowed_opt->FindObject(argv[i]);
        arg="";
        i++;
        while((i<argc) && (!IsOption(allowed_opt,argv[i]))){
          arg+=argv[i];
          arg+=" ";
          i++;
        }
        //Adding to map
        //Checking if unique
        if(!optmap->FindObject(option->GetName())){
          arg.ReplaceAll(" ","");
          optmap->Add(new TObjString(*option), new TObjString(arg));
        } else{
          cout<<"! Error, you specified one option "
              <<option->GetName()<<" 2 times"<<endl;
          abort();
        }
        if(i<argc) i--;
      }
    }
  }
  delete allowed_opt;
  return 0;
}

Bool_t IsOption(TList *optlist, TString word){
  if(!optlist->FindObject(word)) return kFALSE;
  return kTRUE;
}

//Printing user options
void PrintOptions(TMap* optmap){
  cout<<"Application will be executed with next options: "<<endl;
  TMapIter *n= (TMapIter *)optmap->MakeIterator();
  TObject *o; while( (o=(TObject *) n->Next())) {
    cout<<" option: "<<setw(16)
        <<((TObjString*)o)->GetName()
        <<"    value: "
        <<((TObjString*)optmap->GetValue(o))->GetName()<<endl;      
  } 
}


TH2D* MakeLTmaskHisto(TObject* path_to_rig, Int_t nmasks, TimeVec* tm_vec, Int_t bun){


  //rig bins
  TH1D *h_templ_r = Pa::GetBinning( ((TObjString*)path_to_rig)->GetName() );
  const Int_t nrb = h_templ_r->GetNbinsX();
  Double_t *r_bins = new Double_t[nrb+1];
  for(Int_t k = 0; k<nrb+1; k++) r_bins[k]=h_templ_r->GetBinLowEdge(k+1);
  delete h_templ_r;

  //mask bins
  Double_t *m_bins = new Double_t[nmasks+1];
    for(Int_t k = 0; k<nmasks+1; k++) m_bins[k]=0.5+k;
   
  TString name = "LT_hist_bunch_";
  name += bun;
  name.ReplaceAll(" ", "");
  
  TString title = "Bunch_";
  title += bun;
  title.ReplaceAll(" ","");
  title += ", From ";
  UInt_t start, stop;
  tm_vec->GetStartStop( bun, start, stop );
  title += Pa::PrintTime( start );
  title += " To ";
  title += Pa::PrintTime( stop );

  return new TH2D(name, title, nrb, r_bins, nmasks, m_bins);

}

TimeVec * MakeTimeVec( TObject* path_to_textfile){

  if( !path_to_textfile ){
    cout<<" No path to TimeVec textfile specified. Check --bunch_file option "<<endl;
    abort();
  }

  TString path = ((TObjString*)path_to_textfile)->GetName();

  ifstream bun_list;
  bun_list.open(path.Data());
  if(!bun_list){
    cout<<"! Error, Text file with start-stop abstime: "<<path
        <<"  can't be read. Check --bunch_file option"<<endl;
    abort();
  }
  TString start_s,stop_s;


  TimeVec *v_out = new TimeVec();

  while (!bun_list.eof()) {
    bun_list>>start_s;
    bun_list>>stop_s;
    if(start_s != ""){
      v_out->Add( start_s.Atoll(), stop_s.Atoll());
    }
  }

  return v_out;
}


TList * DoMasksProfile( TObject* inp_path_o, TObject* inp_pattern, TObject* path_to_rig, TimeVec* tm_vec, Double_t cut_f, Double_t b_min, Double_t b_max ){

  if( !inp_path_o ){
    cout<<"No input directory defined, check --inp_dir option"<<endl;
    abort();
  }
  if( !inp_pattern ){
    cout<<"No input pattern defined, check --inp_pattern option"<<endl;
    abort();
  }

  if( !path_to_rig ){
    cout<<" No path to Rigidity bins textfile specified. Check --rigbin_file option "<<endl;
    abort();
  }
  
  //template for LT cutoff profile
  TH1D* lt_templ = Pa::GetBinning( ((TObjString*)path_to_rig)->GetName() );

  
  //Here I'm making chain and define some variables
  Double_t babs, svlcutoff, ltime;
  UInt_t abstime;
  UInt_t mask[12];
 
  TChain* ch = new TChain("mask_vk_tree");

  ch->SetBranchAddress("abstime",&abstime);
  ch->SetBranchAddress("babs",&babs);
  ch->SetBranchAddress("ltime",&ltime);
  ch->SetBranchAddress("svlcutoff",&svlcutoff);
  ch->SetBranchAddress("VKmask",&mask);


  //Adding files to chain
  TSystemDirectory* dir = new TSystemDirectory("dir", ((TObjString*)inp_path_o)->GetName());
  TList* ldir = dir->GetListOfFiles();
  TString pattern = ((TObjString*)inp_pattern)->GetName();

  TRegexp re ( pattern, pattern.MaybeWildcard() );
  for(int i=0; i<ldir->GetEntries(); i++) {
      TSystemFile* sysf = (TSystemFile*) ldir->At(i);
      if(sysf->IsDirectory()) continue;
      TString nm = sysf->GetName();
      TString path = sysf->GetTitle();
      TString s = nm;
      if ((nm!=pattern) && s.Index(re) == kNPOS) continue;
      path += sysf->GetName();
      cout<<path<<endl;
      ch->AddFile( path );
  }
  delete ldir;
  delete dir;


  ch->Add( ((TObjString*)inp_path_o)->GetName() );
  cout<<"entries: "<<ch->GetEntries()<<endl;

  //Making general lists
  TList *lout = new TList();
  for(Int_t bun = 0; bun<tm_vec->GetSize(); bun++){
    //fi
    TList *lbun = new TList();
    //first list is list to be saved: TH2D + List of masks
    lbun->Add(new TList() );
    //second list is aux list, not to be saved. Here I'll put TH1Ds for specific mask as function of rigidity cutoff
    lbun->Add(new TList() );
    lout->Add( lbun );
  }

  VKmask *probe_mask = new VKmask();

  for(UInt_t ev = 0; ev < ch->GetEntries(); ev++){
    ch->GetEntry(ev);
    if(ev%1000000==0) cout<<ev<<endl;

    Int_t nb = tm_vec->GetIntervalId( abstime );
    if (nb < 0 ) continue;
    if ( ( babs < b_min) || ( babs > b_max) ) continue; 

    //Filling up current mask
    for(Int_t j =0; j<12; j++) probe_mask->SetMaskRow(j, mask[j]);
    //Compare with already stored masks
    TList *lb = (TList*)lout->At( nb );
    TList *lb_msk =(TList*)lb->At(0);
    TList *lb_hst =(TList*)lb->At(1);
    Int_t id = FindMaskIDinList( lb_msk, probe_mask);
    //if nothing found...
    if(id < 0 ){
      cout<<"-->A new mask... for interval:"<<nb<<endl;
      VKmask *nm = new VKmask();
      for(Int_t j = 0; j<12; j++) nm->SetMaskRow( j, probe_mask->GetMaskRow(j) );
      lb_msk->Add( nm );

      Int_t msk_count = lb_msk->GetEntries();
      TString nm_prof = "bunch_";
      nm_prof += nb;
      nm_prof += "_mask_";
      nm_prof += msk_count;
      TH1D* lt_prof = (TH1D*)lt_templ->Clone( nm_prof );
      Int_t r_b = lt_templ->FindBin( svlcutoff * cut_f );
      for(Int_t i = r_b + 1; i<lt_templ->GetNbinsX()+1; i++)
        lt_prof->AddBinContent( i, ltime );
      lb_hst->Add( lt_prof );   
      cout<<"<--A new mask..."<<endl;
    } else { //if this mask was already there, increment LT
      TH1D* lt_prof = (TH1D*)lb_hst->At( id );
      Int_t r_b = lt_templ->FindBin( svlcutoff * cut_f );
      for(Int_t i = r_b + 1; i<lt_templ->GetNbinsX()+1; i++)
        lt_prof->AddBinContent( i, ltime );
    }

  }

  //After everything we will create TH2D profiles for each bunch
  //TH2D will be added to 1st list (where we have masks) as last object
  for(Int_t bun = 0; bun<tm_vec->GetSize(); bun++){
    TList *lb = (TList*)lout->At(bun);
    TList *lb_msk =(TList*)lb->At(0);
    TList *lb_hst =(TList*)lb->At(1);
    
    const Int_t nmasks = lb_hst->GetEntries(); 
    cout<<"Number of mask for bunch "<<bun<<" is:"<<nmasks<<endl;
    TH2D* h2_mask_prof = MakeLTmaskHisto(path_to_rig, nmasks, tm_vec,  bun);
    for(Int_t msk = 0; msk < nmasks; msk ++ ){
      TH1D* h_prof = (TH1D*)lb_hst->At(msk);
      for(Int_t r_b = 0; r_b<h_prof->GetNbinsX()+1; r_b ++ )
        h2_mask_prof->SetBinContent( r_b+1, msk+1, h_prof->GetBinContent(r_b+1));
    }

    lb_msk->Add( h2_mask_prof );
  }


  return lout;
}

Int_t FindMaskIDinList(TList *lmsk, VKmask* msk){
  
  for(Int_t j = 0; j<lmsk->GetEntries(); j++){
    VKmask *tms = (VKmask*)lmsk->At(j);
    if (*tms==*msk) return j;
  }
  return -1;
}

void SaveResult( TObject* out_path_o, TList* out_list){

 if( !out_path_o ){
    cout<<"No input directory defined, check --out_dir option"<<endl;
    abort();
  }

 TString path = ((TObjString*)out_path_o)->GetName();

 for(Int_t bun = 0; bun < out_list->GetEntries(); bun++){
   TString file = path;
   file += "/out_bun_";
   file += bun;
   file += ".root";
   TList *lb = (TList*)out_list->At(bun);
   TList *lb_msk =(TList*)lb->At(0);
   TFile *fout = new TFile (file, "recreate" );
  
   lb_msk->Write("masks", TObject::kSingleKey);
   fout->Close();
 }

}
1	formato	1.1	#include <iostream>
2			#include <fstream>
3
4			#include <TBenchmark.h>
5			#include <TMap.h>
6			#include <TList.h>
7			#include <TString.h>
8			#include <TObjString.h>
9			#include <TH1D.h>
10			#include <TH2D.h>
11			#include <TSystemDirectory.h>
12			#include <TSystemFile.h>
13			#include <TRegexp.h>
14
15			#include <PamCut/PamCutEnv.h>
16			#include <VKmask.h>
17			#include <TimeVec.h>
18
19
20			//Declaration of all functions
21			Int_t HandleInputPar(Int_t argc, char *argv, TMap optmap);
22			Bool_t IsOption(TList* optlist, TString word);
23			// command line parser
24			void PrintOptions(TMap* optmap);
25			//Making TH2D with LT profiles
26			TH2D* MakeLTmaskHisto(TObject* path_to_rig, Int_t nmasks, TimeVec* tm_vec, Int_t bun);
27			//Make TimeVector from textfile
28			TimeVec * MakeTimeVec( TObject* path_to_textfile);
29			Int_t FindMaskIDinList(TList lmsk, VKmask msk);
30
31
32			//general loop. Returns List of lists to store in root-files
33			TList * DoMasksProfile( TObject* inp_path_o, TObject* inp_pattern, TObject* path_to_rig, TimeVec* tm_vec, Double_t cut_f, Double_t b_min, Double_t b_max );
34
35			void SaveResult( TObject* out_path_o, TList* out_list);
36
37			Int_t main(Int_t argc, char **argv){
38			// Arguments are following:
39			//
40			// --inp_dir /// input directory with root files
41			// --rigbin_file /// file with left edges of rigidity bins
42			// --cutof_factor /// factor to cut in SVL. Typycal value is 1.3
43			// --babs_min /// minimum value of magnetic field. 0.26 to exclude SAA
44			// --babs_max /// maximum value of m.f. usually 1000000 to include everything
45			// --bunch_file /// textfile with bunches with start-stop values
46
47			//Default values
48			Double_t cut_f = 1.3;
49			Double_t babs_min = 0.26;
50			Double_t babs_max = 10000000.;
51
52
53			TMap* optmap = new TMap(); //map with options
54			HandleInputPar(argc,argv, optmap);
55			PrintOptions(optmap);
56
57
58			TBenchmark* bench = new TBenchmark();
59			bench->Start("client");
60
61
62			//Making TimeVec
63
64			TimeVec * tm_vec = MakeTimeVec( optmap->GetValue("--bunch_file") );
65			tm_vec->Print();
66
67
68			TList *lout = DoMasksProfile( optmap->GetValue("--inp_dir"),
69			optmap->GetValue("--inp_pattern"),
70			optmap->GetValue("--rigbin_file"),
71			tm_vec, cut_f, babs_min, babs_max );
72
73
74			SaveResult( optmap->GetValue("--out_dir"), lout);
75
76			// report a benchmark message to the server
77			bench->Stop("client");
78			TString benchresult = " job";
79			benchresult.Remove(0, 1);
80			benchresult += " finished (";
81			TDatime d;
82			benchresult += d.AsSQLString();
83			benchresult += " \| cpu time = ";
84			benchresult += floor((Double_t) bench->GetCpuTime("client"));
85			benchresult += " s = ";
86			benchresult += floor(((Double_t)bench->GetCpuTime("client"))/((Double_t) bench->GetRealTime("client"))*100.);
87			benchresult += " % real time)";
88			cout<<benchresult.Data()<<endl;
89			delete bench;
90
91
92			return 0;
93			}
94
95
96			// Command line parser
97
98			Int_t HandleInputPar(Int_t argc, char *argv, TMap optmap){
99
100
101			TList* allowed_opt = new TList();
102			TString optname = "--inp_dir"; // path with pattern to input files to make chain
103			allowed_opt->Add(new TObjString(optname));
104			optname = "--inp_pattern"; // pattern of files to make chain
105			allowed_opt->Add(new TObjString(optname));
106			optname="--out_dir"; //output path to store root-files
107			allowed_opt->Add(new TObjString(optname));
108			optname = "--rigbin_file"; // full path to textfile with rigidity bins
109			allowed_opt->Add(new TObjString(optname));
110			optname = "--cutoff_factor"; // factor for SVL. Typ. is 1.3
111			allowed_opt->Add(new TObjString(optname));
112			optname = "--babs_min"; // cut to remove SAA and polar regions. Typ 0.26
113			allowed_opt->Add(new TObjString(optname));
114			optname = "--babs_max"; //
115			allowed_opt->Add(new TObjString(optname));
116			optname = "--bunch_file"; // file with start-stop abstime
117			allowed_opt->Add(new TObjString(optname));
118
119
120			TString arg="";
121
122
123
124			if(argc>1){
125			if(!strcmp(argv[1], "-h") \|\| !strcmp(argv[1], "--help") ){
126			return(1);
127			};
128
129			for (Int_t i = 1; i < argc; i++){
130			if(IsOption(allowed_opt,argv[i])){
131			TObjString* option =(TObjString*)allowed_opt->FindObject(argv[i]);
132			arg="";
133			i++;
134			while((i<argc) && (!IsOption(allowed_opt,argv[i]))){
135			arg+=argv[i];
136			arg+=" ";
137			i++;
138			}
139			//Adding to map
140			//Checking if unique
141			if(!optmap->FindObject(option->GetName())){
142			arg.ReplaceAll(" ","");
143			optmap->Add(new TObjString(*option), new TObjString(arg));
144			} else{
145			cout<<"! Error, you specified one option "
146			<<option->GetName()<<" 2 times"<<endl;
147			abort();
148			}
149			if(i<argc) i--;
150			}
151			}
152			}
153			delete allowed_opt;
154			return 0;
155			}
156
157			Bool_t IsOption(TList *optlist, TString word){
158			if(!optlist->FindObject(word)) return kFALSE;
159			return kTRUE;
160			}
161
162			//Printing user options
163			void PrintOptions(TMap* optmap){
164			cout<<"Application will be executed with next options: "<<endl;
165			TMapIter n= (TMapIter )optmap->MakeIterator();
166			TObject o; while( (o=(TObject ) n->Next())) {
167			cout<<" option: "<<setw(16)
168			<<((TObjString*)o)->GetName()
169			<<" value: "
170			<<((TObjString*)optmap->GetValue(o))->GetName()<<endl;
171			}
172			}
173
174
175			TH2D* MakeLTmaskHisto(TObject* path_to_rig, Int_t nmasks, TimeVec* tm_vec, Int_t bun){
176
177
178
179			//rig bins
180			TH1D h_templ_r = Pa::GetBinning( ((TObjString)path_to_rig)->GetName() );
181			const Int_t nrb = h_templ_r->GetNbinsX();
182			Double_t *r_bins = new Double_t[nrb+1];
183			for(Int_t k = 0; k<nrb+1; k++) r_bins[k]=h_templ_r->GetBinLowEdge(k+1);
184			delete h_templ_r;
185
186			//mask bins
187			Double_t *m_bins = new Double_t[nmasks+1];
188			for(Int_t k = 0; k<nmasks+1; k++) m_bins[k]=0.5+k;
189
190			TString name = "LT_hist_bunch_";
191			name += bun;
192			name.ReplaceAll(" ", "");
193
194			TString title = "Bunch_";
195			title += bun;
196			title.ReplaceAll(" ","");
197			title += ", From ";
198			UInt_t start, stop;
199			tm_vec->GetStartStop( bun, start, stop );
200			title += Pa::PrintTime( start );
201			title += " To ";
202			title += Pa::PrintTime( stop );
203
204			return new TH2D(name, title, nrb, r_bins, nmasks, m_bins);
205
206			}
207
208			TimeVec * MakeTimeVec( TObject* path_to_textfile){
209
210			if( !path_to_textfile ){
211			cout<<" No path to TimeVec textfile specified. Check --bunch_file option "<<endl;
212			abort();
213			}
214
215			TString path = ((TObjString*)path_to_textfile)->GetName();
216
217			ifstream bun_list;
218			bun_list.open(path.Data());
219			if(!bun_list){
220			cout<<"! Error, Text file with start-stop abstime: "<<path
221			<<" can't be read. Check --bunch_file option"<<endl;
222			abort();
223			}
224			TString start_s,stop_s;
225
226
227			TimeVec *v_out = new TimeVec();
228
229			while (!bun_list.eof()) {
230			bun_list>>start_s;
231			bun_list>>stop_s;
232			if(start_s != ""){
233			v_out->Add( start_s.Atoll(), stop_s.Atoll());
234			}
235			}
236
237			return v_out;
238			}
239
240
241			TList * DoMasksProfile( TObject* inp_path_o, TObject* inp_pattern, TObject* path_to_rig, TimeVec* tm_vec, Double_t cut_f, Double_t b_min, Double_t b_max ){
242
243			if( !inp_path_o ){
244			cout<<"No input directory defined, check --inp_dir option"<<endl;
245			abort();
246			}
247			if( !inp_pattern ){
248			cout<<"No input pattern defined, check --inp_pattern option"<<endl;
249			abort();
250			}
251
252			if( !path_to_rig ){
253			cout<<" No path to Rigidity bins textfile specified. Check --rigbin_file option "<<endl;
254			abort();
255			}
256
257			//template for LT cutoff profile
258			TH1D* lt_templ = Pa::GetBinning( ((TObjString*)path_to_rig)->GetName() );
259
260
261			//Here I'm making chain and define some variables
262			Double_t babs, svlcutoff, ltime;
263			UInt_t abstime;
264			UInt_t mask[12];
265
266			TChain* ch = new TChain("mask_vk_tree");
267
268			ch->SetBranchAddress("abstime",&abstime);
269			ch->SetBranchAddress("babs",&babs);
270			ch->SetBranchAddress("ltime",&ltime);
271			ch->SetBranchAddress("svlcutoff",&svlcutoff);
272			ch->SetBranchAddress("VKmask",&mask);
273
274
275			//Adding files to chain
276			TSystemDirectory* dir = new TSystemDirectory("dir", ((TObjString*)inp_path_o)->GetName());
277			TList* ldir = dir->GetListOfFiles();
278			TString pattern = ((TObjString*)inp_pattern)->GetName();
279
280			TRegexp re ( pattern, pattern.MaybeWildcard() );
281			for(int i=0; i<ldir->GetEntries(); i++) {
282			TSystemFile* sysf = (TSystemFile*) ldir->At(i);
283			if(sysf->IsDirectory()) continue;
284			TString nm = sysf->GetName();
285			TString path = sysf->GetTitle();
286			TString s = nm;
287			if ((nm!=pattern) && s.Index(re) == kNPOS) continue;
288			path += sysf->GetName();
289			cout<<path<<endl;
290			ch->AddFile( path );
291			}
292			delete ldir;
293			delete dir;
294
295
296
297			ch->Add( ((TObjString*)inp_path_o)->GetName() );
298			cout<<"entries: "<<ch->GetEntries()<<endl;
299
300			//Making general lists
301			TList *lout = new TList();
302			for(Int_t bun = 0; bun<tm_vec->GetSize(); bun++){
303			//fi
304			TList *lbun = new TList();
305			//first list is list to be saved: TH2D + List of masks
306			lbun->Add(new TList() );
307			//second list is aux list, not to be saved. Here I'll put TH1Ds for specific mask as function of rigidity cutoff
308			lbun->Add(new TList() );
309			lout->Add( lbun );
310			}
311
312			VKmask *probe_mask = new VKmask();
313
314			for(UInt_t ev = 0; ev < ch->GetEntries(); ev++){
315			ch->GetEntry(ev);
316			if(ev%1000000==0) cout<<ev<<endl;
317
318			Int_t nb = tm_vec->GetIntervalId( abstime );
319			if (nb < 0 ) continue;
320			if ( ( babs < b_min) \|\| ( babs > b_max) ) continue;
321
322			//Filling up current mask
323			for(Int_t j =0; j<12; j++) probe_mask->SetMaskRow(j, mask[j]);
324			//Compare with already stored masks
325			TList lb = (TList)lout->At( nb );
326			TList lb_msk =(TList)lb->At(0);
327			TList lb_hst =(TList)lb->At(1);
328			Int_t id = FindMaskIDinList( lb_msk, probe_mask);
329			//if nothing found...
330			if(id < 0 ){
331			cout<<"-->A new mask... for interval:"<<nb<<endl;
332			VKmask *nm = new VKmask();
333			for(Int_t j = 0; j<12; j++) nm->SetMaskRow( j, probe_mask->GetMaskRow(j) );
334			lb_msk->Add( nm );
335
336			Int_t msk_count = lb_msk->GetEntries();
337			TString nm_prof = "bunch_";
338			nm_prof += nb;
339			nm_prof += "_mask_";
340			nm_prof += msk_count;
341			TH1D* lt_prof = (TH1D*)lt_templ->Clone( nm_prof );
342			Int_t r_b = lt_templ->FindBin( svlcutoff * cut_f );
343			for(Int_t i = r_b + 1; i<lt_templ->GetNbinsX()+1; i++)
344			lt_prof->AddBinContent( i, ltime );
345			lb_hst->Add( lt_prof );
346			cout<<"<--A new mask..."<<endl;
347			} else { //if this mask was already there, increment LT
348			TH1D* lt_prof = (TH1D*)lb_hst->At( id );
349			Int_t r_b = lt_templ->FindBin( svlcutoff * cut_f );
350			for(Int_t i = r_b + 1; i<lt_templ->GetNbinsX()+1; i++)
351			lt_prof->AddBinContent( i, ltime );
352			}
353
354			}
355
356			//After everything we will create TH2D profiles for each bunch
357			//TH2D will be added to 1st list (where we have masks) as last object
358			for(Int_t bun = 0; bun<tm_vec->GetSize(); bun++){
359			TList lb = (TList)lout->At(bun);
360			TList lb_msk =(TList)lb->At(0);
361			TList lb_hst =(TList)lb->At(1);
362
363			const Int_t nmasks = lb_hst->GetEntries();
364			cout<<"Number of mask for bunch "<<bun<<" is:"<<nmasks<<endl;
365			TH2D* h2_mask_prof = MakeLTmaskHisto(path_to_rig, nmasks, tm_vec, bun);
366			for(Int_t msk = 0; msk < nmasks; msk ++ ){
367			TH1D* h_prof = (TH1D*)lb_hst->At(msk);
368			for(Int_t r_b = 0; r_b<h_prof->GetNbinsX()+1; r_b ++ )
369			h2_mask_prof->SetBinContent( r_b+1, msk+1, h_prof->GetBinContent(r_b+1));
370			}
371
372			lb_msk->Add( h2_mask_prof );
373			}
374
375
376			return lout;
377			}
378
379			Int_t FindMaskIDinList(TList lmsk, VKmask msk){
380
381			for(Int_t j = 0; j<lmsk->GetEntries(); j++){
382			VKmask tms = (VKmask)lmsk->At(j);
383			if (tms==msk) return j;
384			}
385			return -1;
386			}
387
388			void SaveResult( TObject* out_path_o, TList* out_list){
389
390			if( !out_path_o ){
391			cout<<"No input directory defined, check --out_dir option"<<endl;
392			abort();
393			}
394
395			TString path = ((TObjString*)out_path_o)->GetName();
396
397			for(Int_t bun = 0; bun < out_list->GetEntries(); bun++){
398			TString file = path;
399			file += "/out_bun_";
400			file += bun;
401			file += ".root";
402			TList lb = (TList)out_list->At(bun);
403			TList lb_msk =(TList)lb->At(0);
404			TFile *fout = new TFile (file, "recreate" );
405
406			lb_msk->Write("masks", TObject::kSingleKey);
407			fout->Close();
408			}
409
410			}