| 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",<ime); |
| 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 |
} |
Parent Directory
| 
Revision Log