18 #include "TObjString.h"
20 #include "TObjArray.h"
21 #include "TObjString.h"
143 unsigned int nNull = 0;
145 bool ok = (
size > 0);
151 bool ok = (
size > 0);
157 cout << i1 <<
"-" << i2 << endl;
179 std::vector<unsigned int> hashV, indexV;
180 unsigned int nTot = 0;
184 if (!history)
continue;
185 for (
unsigned int j = 0; j < history->
nData(); j++) {
188 if (
data.isDisconnected())
continue;
189 if (nTot % 10000 == 0) cout << nTot << endl;
190 if (
data.energy() > eCut) {
191 cout <<
"E = " <<
data.energy() <<
" " <<
i <<
" " << j << endl;
201 for (
unsigned int i = 0;
i < hashV.size();
i++) {
206 cout << hashV.size() <<
"/" << nTot << endl;
215 if (!history)
continue;
217 cout <<
"Invalid LArSamplesHistory at hash = " <<
i << endl;
228 std::vector<const Interface*> interfaces;
229 interfaces.push_back(
this);
230 interfaces.push_back(&
other);
238 for (
unsigned int i = 0;
i < interfaces.size();
i++)
239 accessors.push_back(&interfaces[
i]->accessor());
246 std::vector<const Interface*> interfaces;
247 interfaces.push_back(
this);
248 interfaces.push_back(&
other);
256 for (
unsigned int i = 0;
i < interfaces.size();
i++)
257 accessors.push_back(&interfaces[
i]->accessor());
266 if (!
multi)
return nullptr;
267 std::vector<const Interface*> justOne;
268 justOne.push_back(
multi);
275 if (!
multi)
return nullptr;
276 std::vector<const Interface*> justOne;
277 justOne.push_back(
multi);
286 TObjArray*
list = filters.Tokenize(
",;");
287 if (
list->GetEntries() == 0) {
288 cout <<
"No filtering specified, exiting.";
293 for (
int k = 0;
k <
list->GetEntries();
k++) {
294 TObjString* tobs = (TObjString*)(
list->At(
k));
295 TObjArray*
items = tobs->String().Tokenize(
":");
296 if (
items->GetEntries() != 2) {
297 cout <<
"Invalid filter entry " << tobs->String() <<
", exiting." << endl;
302 TString
params = ((TObjString*)(
items->At(0)))->String();
303 TString
suffix = ((TObjString*)(
items->At(1)))->String();
306 cout <<
"---" << endl;
314 if (!tweak.
set(tweaks))
return 0;
317 if (!
multi)
return 0;
323 std::vector<MultiTreeAccessor*> filtered_mts = mt->
filterComponents(filterList, tweak);
324 if (filtered_mts.size() != filterList.
size()){
329 cout <<
"Component filtering done!" << endl;
333 for (
unsigned int f = 0;
f < filtered_mts.size();
f++) {
334 std::vector<TString>
files;
335 for (
unsigned int i = 0;
i < filtered_mts[
f]->nAccessors();
i++) {
337 cout <<
"Added " <<
files.back() << endl;
339 delete filtered_mts[
f];
342 std::vector<const Interface*> justOne;
343 justOne.push_back(filtered_multi);
355 if (!
f.set(
sel))
return nullptr;
358 if (!tweak.
set(tweaks))
return nullptr;
363 if (not pAccess)
return nullptr;
365 if (newFN !=
"") thisFN = std::move(newFN);
367 return filter(
f, tweak, thisFN);
376 return rootName +
"_" +
suffix +
".root";
400 return filter(
f, tw, newFileName);
413 if (
info->region() != region)
continue;
427 if (
info->feb() != feb)
continue;
442 if (
ft >= 0 &&
info->feedThrough() !=
ft)
continue;
443 if (slot >= 0 &&
info->slot() != slot)
continue;
456 if (!
f.set(
sel))
return nullptr;
458 std::vector<TString> vars;
459 std::vector<DataFuncSet> funcs;
460 std::vector<DataFuncArgs>
args;
462 cout <<
"Invalid variable specification " <<
var << endl;
466 TString
title = vars[0];
469 TH1D*
h =
m.dist(funcs[0],
args[0], vars[0],
nBins, xMin, xMax,
470 title, vars[0],
"digits",
f);
471 if (!
h)
return nullptr;
478 int nBinsY,
double yMin,
double yMax,
479 const TString&
sel,
const TString&
opt)
const
483 if (!
f.set(
sel))
return nullptr;
485 std::vector<TString> vars;
486 std::vector<DataFuncSet> funcs;
487 std::vector<DataFuncArgs>
args;
489 cout <<
"Invalid variable specification " <<
varList << endl;
493 TString
title = vars[1] +
" vs. " + vars[0];
496 TH2D*
h =
m.dist(funcs[0],
args[0], funcs[1],
args[1], vars[0] +
"_" + vars[1],
497 nBinsX, xMin, xMax, nBinsY, yMin, yMax,
498 title, vars[0], vars[1],
f);
499 if (!
h)
return nullptr;
506 const TString&
sel,
const TString&
opt,
511 if (!
f.set(
sel))
return nullptr;
513 std::vector<TString> vars;
514 std::vector<DataFuncSet> funcs;
515 std::vector<DataFuncArgs>
args;
517 cout <<
"Invalid variable specification " <<
var << endl;
526 if (!
h)
return nullptr;
533 const TString&
sel,
const TString&
opt,
538 if (!
f.set(
sel))
return nullptr;
540 std::vector<TString> vars;
541 std::vector<DataFuncSet> funcs;
542 std::vector<DataFuncArgs>
args;
544 cout <<
"Invalid variable specification " <<
var << endl;
552 if (!
h)
return nullptr;
562 if (!
f.set(
sel))
return 0;
571 if (!
f.set(
sel))
return 0;
579 if (!history)
return false;
587 if (!
f.set(
sel))
return false;
590 if (!history)
continue;
594 if (hDesc ==
"")
continue;
595 cout << Form(
"Hash = %-5d : ",
i) << hDesc
596 <<
"-----------------------------------------------------------------------------"
606 if (!
f.set(
sel))
return false;
607 std::map< std::pair<unsigned int, unsigned int>,
unsigned int > eventCells;
608 std::map< std::pair<unsigned int, unsigned int>,
double > eventEnergy;
612 if ((
i+1) % 50000 == 0) cout <<
"Cell # " <<
i+1 <<
"/" <<
nChannels() << endl;
614 if (!history)
continue;
615 for (
unsigned int j = 0; j < history->
nData(); j++) {
616 std::pair<unsigned int, unsigned int>
ev = std::make_pair(history->
data(j)->
run(), history->
data(j)->
event());
622 for (
unsigned int i = 0;
i <
nEvents();
i++) {
624 std::pair<unsigned int, unsigned int>
id(evtData->
run(), evtData->
event());
627 printout += Form(
" : %6d LAr hits, %7.1f MeV", eventCells[
id], eventEnergy[
id]);
628 cout << printout << endl;
637 std::map<unsigned int, unsigned int>
events;
644 for (
unsigned int i = 0;
i <
nRuns();
i++) {
648 printout += Form(
" : %6d events",
events[rData->
run()]);
649 cout << printout << endl;
660 if (!
f.set(
sel))
return 0;
662 std::vector<TString> vars;
663 std::vector<DataFuncSet> funcs;
664 std::vector<DataFuncArgs>
args;
666 cout <<
"Invalid variable specification " <<
varList << endl;
670 TVectorD
mean(vars.size()), meanErr(vars.size());
671 TMatrixD
covMatrix(vars.size(), vars.size()), covMatrixErr(vars.size(), vars.size());
675 cout <<
"---------------------------" << endl;
676 for (
unsigned int i = 0;
i < vars.size();
i++)
677 cout << Form(
"| %10s | %-9.4g |", vars[
i].
Data(),
mean(
i)) << endl;
678 cout <<
"---------------------------" << endl << endl;
681 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
" |";
682 cout << endl <<
"--------------";
683 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"-------------";
685 for (
unsigned int i1 = 0; i1 < vars.size(); i1++) {
686 cout << Form(
"| %10s |", vars[i1].
Data());
687 for (
unsigned int i2 = 0; i2 < vars.size(); i2++)
688 cout << Form(
" %-9.4g |",
covMatrix(i1, i2));
691 cout <<
"--------------";
692 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"-------------";
696 cout <<
"---------------------------" << endl;
697 for (
unsigned int i = 0;
i < vars.size();
i++)
698 cout << Form(
"| %10s | %-9.4g +/- %-9.4g |", vars[
i].
Data(),
mean(
i), meanErr(
i)) << endl;
699 cout <<
"---------------------------" << endl << endl;
702 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
" |";
703 cout << endl <<
"--------------";
704 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"---------------------------";
706 for (
unsigned int i1 = 0; i1 < vars.size(); i1++) {
707 cout << Form(
"| %10s |", vars[i1].
Data());
708 for (
unsigned int i2 = 0; i2 < vars.size(); i2++)
709 cout << Form(
" %-9.4g +/- %-9.4g |",
covMatrix(i1, i2), covMatrixErr(i1, i2));
712 cout <<
"--------------";
713 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"---------------------------";
725 if (!otherCell)
continue;
726 if (
cell.position().DeltaR(otherCell->
position()) > dRCut) {
delete otherCell;
continue; }
738 if (!
cell)
return true;
740 if (
layer < 0)
return true;
741 std::vector<unsigned int> allHashes;
749 if (
layer == -2) {
hashes = allHashes;
return true; }
750 for (
unsigned int h : allHashes) {
775 if (!history)
continue;
776 const Data* dataForEvent = history->data_for_event(
event);
777 if (dataForEvent)
data.push_back(
new Data(*dataForEvent));
785 std::vector<float*> floatVars;
786 std::vector<int*> intVars;
787 std::vector<std::vector<float>*> floatVects;
788 std::vector<std::vector<int>*> intVects;
789 std::map<TString, unsigned int> varIndex;
791 std::vector<TString> vars;
792 std::vector<DataFuncSet> funcs;
793 std::vector<DataFuncArgs>
args;
796 cout <<
"Making trees..." << endl;
798 TFile* flatFile = TFile::Open(
fileName +
"_tmpFlatFile.root",
"RECREATE");
799 TTree* flatTree =
new TTree(
"flatTree",
"Flat tree");
801 TFile* eventFile = TFile::Open(
fileName,
"RECREATE");
802 TTree* eventTree =
new TTree(
"eventTree",
"Event tree");
804 for (
unsigned int j = 0; j < vars.size(); j++) {
805 unsigned int index = 0;
806 if (funcs[j].isNull())
return false;
807 if (funcs[j].isInt()) {
808 int* varCont =
new int(0);
809 std::vector<int>* vectCont =
new std::vector<int>();
810 index = intVars.size();
811 intVars.push_back(varCont);
812 intVects.push_back(vectCont);
813 flatTree->Branch(vars[j], varCont);
814 eventTree->Branch(vars[j], vectCont);
817 float* varCont =
new float(0);
818 std::vector<float>* vectCont =
new std::vector<float>();
819 index = floatVars.size();
820 floatVars.push_back(varCont);
821 floatVects.push_back(vectCont);
822 flatTree->Branch(vars[j], varCont);
823 eventTree->Branch(vars[j], vectCont);
825 varIndex[vars[j]] =
index;
826 cout << vars[j] <<
" -> " <<
index << endl;
828 std::map< unsigned int, std::map< unsigned int, std::vector<long long> > > runEventIndices;
829 cout <<
"Making flat ntuple" << endl;
830 unsigned int count = 0;
834 if (
count % 100 == 0) cout <<
"Processing entry " <<
count <<
" (hash = " << iter.pos() <<
"), size = " <<
hist->nData() << endl;
835 for (
unsigned int k = 0;
k <
hist->nData();
k++) {
836 for (
unsigned int j = 0; j < vars.size(); j++) {
837 if (funcs[j].isInt())
838 *intVars[varIndex[vars[j]]] =
int(funcs[j].intVal(*
hist->data(
k),
args[j]));
840 *floatVars[varIndex[vars[j]]] = funcs[j].doubleVal(*
hist->data(
k),
args[j]);
842 runEventIndices[
hist->data(
k)->run()][
hist->data(
k)->event()].push_back(flatTree->GetEntries());
847 cout <<
"Making event tuple" << endl;
849 for (
const auto&
run : runEventIndices) {
851 cout <<
"Processing run " <<
run.first <<
" (" <<
runCount <<
" of " << runEventIndices.size() <<
")" << endl;
852 unsigned int eventCount = 0;
853 for (
const auto&
event :
run.second) {
855 if (eventCount % 1000 == 0)
856 cout <<
" processing event " <<
event.first <<
" (" << eventCount <<
" of " <<
run.second.size() <<
"), size = " <<
event.second.size() << endl;
857 for (
unsigned int j = 0; j < intVects.size(); j++) intVects[j]->
clear();
858 for (
unsigned int j = 0; j < floatVects.size(); j++) floatVects[j]->
clear();
860 flatTree->GetEntry(
index);
861 for (
unsigned int j = 0; j < vars.size(); j++) {
862 if (funcs[j].isInt())
863 intVects[varIndex[vars[j]]]->push_back(*intVars[varIndex[vars[j]]]);
865 floatVects[varIndex[vars[j]]]->push_back(*floatVars[varIndex[vars[j]]]);
872 cout <<
"Writing data..." << endl;
878 cout <<
"Cleaning up..." << endl;
884 for (
unsigned int j = 0; j < intVects.size(); j++)
delete intVects[j];
885 for (
unsigned int j = 0; j < floatVects.size(); j++)
delete floatVects[j];
886 for (
unsigned int j = 0; j < intVars.size(); j++)
delete intVars[j];
887 for (
unsigned int j = 0; j < floatVars.size(); j++)
delete floatVars[j];
889 cout <<
"Done!" << endl;