20 #include "TObjString.h"
22 #include "TObjArray.h"
23 #include "TObjString.h"
154 unsigned int nNull = 0;
156 bool ok = (
size > 0);
162 bool ok = (
size > 0);
168 cout << i1 <<
"-" << i2 << endl;
190 std::vector<unsigned int> hashV, indexV;
191 unsigned int nTot = 0;
195 if (!history)
continue;
196 for (
unsigned int j = 0; j < history->
nData(); j++) {
199 if (
data.isDisconnected())
continue;
200 if (nTot % 10000 == 0) cout << nTot << endl;
201 if (
data.energy() > eCut) {
202 cout <<
"E = " <<
data.energy() <<
" " <<
i <<
" " << j << endl;
212 for (
unsigned int i = 0;
i < hashV.size();
i++) {
217 cout << hashV.size() <<
"/" << nTot << endl;
226 if (!history)
continue;
228 cout <<
"Invalid LArSamplesHistory at hash = " <<
i << endl;
239 std::vector<const Interface*> interfaces;
240 interfaces.push_back(
this);
241 interfaces.push_back(&
other);
249 for (
unsigned int i = 0;
i < interfaces.size();
i++)
250 accessors.push_back(&interfaces[
i]->accessor());
257 std::vector<const Interface*> interfaces;
258 interfaces.push_back(
this);
259 interfaces.push_back(&
other);
267 for (
unsigned int i = 0;
i < interfaces.size();
i++)
268 accessors.push_back(&interfaces[
i]->accessor());
277 if (!
multi)
return nullptr;
278 std::vector<const Interface*> justOne;
279 justOne.push_back(
multi);
286 if (!
multi)
return nullptr;
287 std::vector<const Interface*> justOne;
288 justOne.push_back(
multi);
297 TObjArray*
list = filters.Tokenize(
",;");
298 if (
list->GetEntries() == 0) {
299 cout <<
"No filtering specified, exiting.";
304 for (
int k = 0;
k <
list->GetEntries();
k++) {
305 TObjString* tobs = (TObjString*)(
list->At(
k));
306 TObjArray*
items = tobs->String().Tokenize(
":");
307 if (
items->GetEntries() != 2) {
308 cout <<
"Invalid filter entry " << tobs->String() <<
", exiting." << endl;
313 TString
params = ((TObjString*)(
items->At(0)))->String();
314 TString
suffix = ((TObjString*)(
items->At(1)))->String();
317 cout <<
"---" << endl;
325 if (!tweak.
set(tweaks))
return 0;
328 if (!
multi)
return 0;
334 std::vector<MultiTreeAccessor*> filtered_mts = mt->
filterComponents(filterList, tweak);
335 if (filtered_mts.size() != filterList.
size()){
340 cout <<
"Component filtering done!" << endl;
344 for (
unsigned int f = 0;
f < filtered_mts.size();
f++) {
345 std::vector<TString>
files;
346 for (
unsigned int i = 0;
i < filtered_mts[
f]->nAccessors();
i++) {
348 cout <<
"Added " <<
files.back() << endl;
350 delete filtered_mts[
f];
353 std::vector<const Interface*> justOne;
354 justOne.push_back(filtered_multi);
366 if (!
f.set(
sel))
return nullptr;
369 if (!tweak.
set(tweaks))
return nullptr;
374 if (not pAccess)
return nullptr;
376 if (newFN !=
"") thisFN = std::move(newFN);
378 return filter(
f, tweak, thisFN);
387 return rootName +
"_" +
suffix +
".root";
411 return filter(
f, tw, newFileName);
424 if (
info->region() != region)
continue;
438 if (
info->feb() != feb)
continue;
453 if (
ft >= 0 &&
info->feedThrough() !=
ft)
continue;
454 if (slot >= 0 &&
info->slot() != slot)
continue;
467 if (!
f.set(
sel))
return nullptr;
469 std::vector<TString> vars;
470 std::vector<DataFuncSet> funcs;
471 std::vector<DataFuncArgs>
args;
473 cout <<
"Invalid variable specification " <<
var << endl;
477 TString
title = vars[0];
480 TH1D*
h =
m.dist(funcs[0],
args[0], vars[0],
nBins, xMin, xMax,
481 title, vars[0],
"digits",
f);
482 if (!
h)
return nullptr;
489 int nBinsY,
double yMin,
double yMax,
490 const TString&
sel,
const TString&
opt)
const
494 if (!
f.set(
sel))
return nullptr;
496 std::vector<TString> vars;
497 std::vector<DataFuncSet> funcs;
498 std::vector<DataFuncArgs>
args;
500 cout <<
"Invalid variable specification " <<
varList << endl;
504 TString
title = vars[1] +
" vs. " + vars[0];
507 TH2D*
h =
m.dist(funcs[0],
args[0], funcs[1],
args[1], vars[0] +
"_" + vars[1],
508 nBinsX, xMin, xMax, nBinsY, yMin, yMax,
509 title, vars[0], vars[1],
f);
510 if (!
h)
return nullptr;
517 const TString&
sel,
const TString&
opt,
522 if (!
f.set(
sel))
return nullptr;
524 std::vector<TString> vars;
525 std::vector<DataFuncSet> funcs;
526 std::vector<DataFuncArgs>
args;
528 cout <<
"Invalid variable specification " <<
var << endl;
537 if (!
h)
return nullptr;
544 const TString&
sel,
const TString&
opt,
549 if (!
f.set(
sel))
return nullptr;
551 std::vector<TString> vars;
552 std::vector<DataFuncSet> funcs;
553 std::vector<DataFuncArgs>
args;
555 cout <<
"Invalid variable specification " <<
var << endl;
563 if (!
h)
return nullptr;
573 if (!
f.set(
sel))
return 0;
582 if (!
f.set(
sel))
return 0;
590 if (!history)
return false;
598 if (!
f.set(
sel))
return false;
601 if (!history)
continue;
605 if (hDesc ==
"")
continue;
606 cout << Form(
"Hash = %-5d : ",
i) << hDesc
607 <<
"-----------------------------------------------------------------------------"
617 if (!
f.set(
sel))
return false;
618 std::map< std::pair<unsigned int, unsigned int>,
unsigned int > eventCells;
619 std::map< std::pair<unsigned int, unsigned int>,
double > eventEnergy;
623 if ((
i+1) % 50000 == 0) cout <<
"Cell # " <<
i+1 <<
"/" <<
nChannels() << endl;
625 if (!history)
continue;
626 for (
unsigned int j = 0; j < history->
nData(); j++) {
627 std::pair<unsigned int, unsigned int>
ev = std::make_pair(history->
data(j)->
run(), history->
data(j)->
event());
633 for (
unsigned int i = 0;
i <
nEvents();
i++) {
635 std::pair<unsigned int, unsigned int>
id(evtData->
run(), evtData->
event());
638 printout += Form(
" : %6d LAr hits, %7.1f MeV", eventCells[
id], eventEnergy[
id]);
639 cout << printout << endl;
648 std::map<unsigned int, unsigned int>
events;
655 for (
unsigned int i = 0;
i <
nRuns();
i++) {
659 printout += Form(
" : %6d events",
events[rData->
run()]);
660 cout << printout << endl;
671 if (!
f.set(
sel))
return 0;
673 std::vector<TString> vars;
674 std::vector<DataFuncSet> funcs;
675 std::vector<DataFuncArgs>
args;
677 cout <<
"Invalid variable specification " <<
varList << endl;
681 TVectorD
mean(vars.size()), meanErr(vars.size());
682 TMatrixD
covMatrix(vars.size(), vars.size()), covMatrixErr(vars.size(), vars.size());
686 cout <<
"---------------------------" << endl;
687 for (
unsigned int i = 0;
i < vars.size();
i++)
688 cout << Form(
"| %10s | %-9.4g |", vars[
i].
Data(),
mean(
i)) << endl;
689 cout <<
"---------------------------" << endl << endl;
692 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
" |";
693 cout << endl <<
"--------------";
694 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"-------------";
696 for (
unsigned int i1 = 0; i1 < vars.size(); i1++) {
697 cout << Form(
"| %10s |", vars[i1].
Data());
698 for (
unsigned int i2 = 0; i2 < vars.size(); i2++)
699 cout << Form(
" %-9.4g |",
covMatrix(i1, i2));
702 cout <<
"--------------";
703 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"-------------";
707 cout <<
"---------------------------" << endl;
708 for (
unsigned int i = 0;
i < vars.size();
i++)
709 cout << Form(
"| %10s | %-9.4g +/- %-9.4g |", vars[
i].
Data(),
mean(
i), meanErr(
i)) << endl;
710 cout <<
"---------------------------" << endl << endl;
713 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
" |";
714 cout << endl <<
"--------------";
715 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"---------------------------";
717 for (
unsigned int i1 = 0; i1 < vars.size(); i1++) {
718 cout << Form(
"| %10s |", vars[i1].
Data());
719 for (
unsigned int i2 = 0; i2 < vars.size(); i2++)
720 cout << Form(
" %-9.4g +/- %-9.4g |",
covMatrix(i1, i2), covMatrixErr(i1, i2));
723 cout <<
"--------------";
724 for (
unsigned int i = 0;
i < vars.size();
i++) cout <<
"---------------------------";
736 if (!otherCell)
continue;
737 if (
cell.position().DeltaR(otherCell->
position()) > dRCut) {
delete otherCell;
continue; }
749 if (!
cell)
return true;
751 if (
layer < 0)
return true;
752 std::vector<unsigned int> allHashes;
760 if (
layer == -2) {
hashes = allHashes;
return true; }
761 for (
unsigned int h : allHashes) {
786 if (!history)
continue;
787 const Data* dataForEvent = history->data_for_event(
event);
788 if (dataForEvent)
data.push_back(
new Data(*dataForEvent));
796 std::vector<float*> floatVars;
797 std::vector<int*> intVars;
798 std::vector<std::vector<float>*> floatVects;
799 std::vector<std::vector<int>*> intVects;
800 std::map<TString, unsigned int> varIndex;
802 std::vector<TString> vars;
803 std::vector<DataFuncSet> funcs;
804 std::vector<DataFuncArgs>
args;
807 cout <<
"Making trees..." << endl;
809 TFile* flatFile = TFile::Open(
fileName +
"_tmpFlatFile.root",
"RECREATE");
810 TTree* flatTree =
new TTree(
"flatTree",
"Flat tree");
812 TFile* eventFile = TFile::Open(
fileName,
"RECREATE");
813 TTree* eventTree =
new TTree(
"eventTree",
"Event tree");
815 for (
unsigned int j = 0; j < vars.size(); j++) {
816 unsigned int index = 0;
817 if (funcs[j].isNull())
return false;
818 if (funcs[j].isInt()) {
819 int* varCont =
new int(0);
820 std::vector<int>* vectCont =
new std::vector<int>();
821 index = intVars.size();
822 intVars.push_back(varCont);
823 intVects.push_back(vectCont);
824 flatTree->Branch(vars[j], varCont);
825 eventTree->Branch(vars[j], vectCont);
828 float* varCont =
new float(0);
829 std::vector<float>* vectCont =
new std::vector<float>();
830 index = floatVars.size();
831 floatVars.push_back(varCont);
832 floatVects.push_back(vectCont);
833 flatTree->Branch(vars[j], varCont);
834 eventTree->Branch(vars[j], vectCont);
836 varIndex[vars[j]] =
index;
837 cout << vars[j] <<
" -> " <<
index << endl;
839 std::map< unsigned int, std::map< unsigned int, std::vector<long long> > > runEventIndices;
840 cout <<
"Making flat ntuple" << endl;
841 unsigned int count = 0;
845 if (
count % 100 == 0) cout <<
"Processing entry " <<
count <<
" (hash = " <<
iter.pos() <<
"), size = " <<
hist->nData() << endl;
846 for (
unsigned int k = 0;
k <
hist->nData();
k++) {
847 for (
unsigned int j = 0; j < vars.size(); j++) {
848 if (funcs[j].isInt())
849 *intVars[varIndex[vars[j]]] =
int(funcs[j].intVal(*
hist->data(
k),
args[j]));
851 *floatVars[varIndex[vars[j]]] = funcs[j].doubleVal(*
hist->data(
k),
args[j]);
853 runEventIndices[
hist->data(
k)->run()][
hist->data(
k)->event()].push_back(flatTree->GetEntries());
858 cout <<
"Making event tuple" << endl;
860 for (
const auto&
run : runEventIndices) {
862 cout <<
"Processing run " <<
run.first <<
" (" <<
runCount <<
" of " << runEventIndices.size() <<
")" << endl;
863 unsigned int eventCount = 0;
864 for (
const auto&
event :
run.second) {
866 if (eventCount % 1000 == 0)
867 cout <<
" processing event " <<
event.first <<
" (" << eventCount <<
" of " <<
run.second.size() <<
"), size = " <<
event.second.size() << endl;
868 for (
unsigned int j = 0; j < intVects.size(); j++) intVects[j]->
clear();
869 for (
unsigned int j = 0; j < floatVects.size(); j++) floatVects[j]->
clear();
871 flatTree->GetEntry(
index);
872 for (
unsigned int j = 0; j < vars.size(); j++) {
873 if (funcs[j].isInt())
874 intVects[varIndex[vars[j]]]->push_back(*intVars[varIndex[vars[j]]]);
876 floatVects[varIndex[vars[j]]]->push_back(*floatVars[varIndex[vars[j]]]);
883 cout <<
"Writing data..." << endl;
889 cout <<
"Cleaning up..." << endl;
895 for (
unsigned int j = 0; j < intVects.size(); j++)
delete intVects[j];
896 for (
unsigned int j = 0; j < floatVects.size(); j++)
delete floatVects[j];
897 for (
unsigned int j = 0; j < intVars.size(); j++)
delete intVars[j];
898 for (
unsigned int j = 0; j < floatVars.size(); j++)
delete floatVars[j];
900 cout <<
"Done!" << endl;