41int main(
int argc,
char* argv[])
43 std::string inputFileName(
"AO2D.root");
44 std::string outputFileName(
"AO2D_thinned.root");
46 bool bOverwrite =
false;
47 int compression = 505;
51 const char*
const short_opts =
"i:o:KOh";
52 static struct option long_options[] = {
53 {
"input", required_argument,
nullptr,
'i'},
54 {
"output", required_argument,
nullptr,
'o'},
55 {
"overwrite", no_argument,
nullptr,
'O'},
56 {
"compression", required_argument,
nullptr,
'c'},
57 {
"help", no_argument,
nullptr,
'h'},
58 {
nullptr, 0,
nullptr, 0}};
61 const auto opt = getopt_long(argc, argv, short_opts, long_options, &option_index);
67 inputFileName = optarg;
70 outputFileName = optarg;
74 printf(
"Overwriting existing output file if existing\n");
77 compression = atoi(optarg);
82 printf(
"AO2D thinning tool. Options: \n");
83 printf(
" --input/-i <inputfile.root> Contains input file path to the file to be thinned. Default: %s\n", inputFileName.c_str());
84 printf(
" --output/-o <outputfile.root> Target output ROOT file. Default: %s\n", outputFileName.c_str());
85 printf(
" --compression/-c <compression id> ROOT compression algorithm / level. Default: %d\n", compression);
87 printf(
" Optional Arguments:\n");
88 printf(
" --overwrite/-O Overwrite existing output file\n");
93 printf(
"AOD reduction started with:\n");
94 printf(
" Input file: %s\n", inputFileName.c_str());
95 printf(
" Ouput file name: %s\n", outputFileName.c_str());
100 auto outputFile = TFile::Open(outputFileName.c_str(), (bOverwrite) ?
"RECREATE" :
"CREATE",
"", compression);
101 if (outputFile ==
nullptr) {
102 printf(
"Error: File %s exists or cannot be created!\n", outputFileName.c_str());
105 TDirectory* outputDir =
nullptr;
107 if (inputFileName.find(
"alien:") == 0 && !gGrid) {
108 printf(
"Connecting to AliEn...");
109 TGrid::Connect(
"alien:");
112 auto inputFile = TFile::Open(inputFileName.c_str());
114 printf(
"Error: Could not open input file %s.\n", inputFileName.c_str());
118 TList* keyList = inputFile->GetListOfKeys();
121 for (
auto key1 : *keyList) {
123 if (((TObjString*)key1)->GetString().EqualTo(
"metaData")) {
124 auto metaData = (TMap*)inputFile->Get(
"metaData");
126 metaData->Write(
"metaData", TObject::kSingleKey);
130 if (((TObjString*)key1)->GetString().EqualTo(
"parentFiles")) {
131 auto parentFiles = (TMap*)inputFile->Get(
"parentFiles");
133 parentFiles->Write(
"parentFiles", TObject::kSingleKey);
137 if (!((TObjString*)key1)->GetString().BeginsWith(
"DF_")) {
141 auto dfName = ((TObjString*)key1)->GetString().Data();
143 printf(
" Processing folder %s\n", dfName);
144 auto folder = (TDirectoryFile*)inputFile->Get(dfName);
145 auto treeList = folder->GetListOfKeys();
149 for (
auto i = 0;
i < treeList->GetEntries(); ++
i) {
150 TKey* ki = (TKey*)treeList->At(
i);
151 for (
int j =
i + 1;
j < treeList->GetEntries(); ++
j) {
152 TKey* kj = (TKey*)treeList->At(
j);
153 if (std::strcmp(ki->GetName(), kj->GetName()) == 0 && std::strcmp(ki->GetTitle(), kj->GetTitle()) == 0) {
154 if (ki->GetCycle() < kj->GetCycle()) {
155 printf(
" *** FATAL *** we had ordered the keys, first cycle should be higher, please check");
159 treeList->Remove(kj);
170 TString v0Name{
"O2v0_???"}, trkExtraName{
"O2trackextra*"};
171 TRegexp v0Re(v0Name, kTRUE), trkExtraRe(trkExtraName, kTRUE);
172 bool hasTrackQA{
false};
173 TString trackQAName{
"O2trackqa*"};
174 TRegexp trackQARe(trackQAName, kTRUE);
175 for (
TObject* obj : *treeList) {
176 TString
st = obj->GetName();
177 if (
st.Index(v0Re) != kNPOS) {
179 }
else if (
st.Index(trkExtraRe) != kNPOS) {
181 }
else if (
st.Index(trackQARe) != kNPOS) {
188 auto v0Entry = (
TObject*)treeList->FindObject(v0Name);
189 treeList->Remove(v0Entry);
190 treeList->AddFirst(v0Entry);
193 auto trackExtraTree = (TTree*)inputFile->Get(Form(
"%s/%s", dfName, trkExtraName.Data()));
194 if (trackExtraTree ==
nullptr) {
195 printf(
"%s table not found\n", trkExtraName.Data());
199 auto track_iu = (TTree*)inputFile->Get(Form(
"%s/%s", dfName,
"O2track_iu"));
200 if (track_iu ==
nullptr) {
201 printf(
"O2track_iu table not found\n");
205 auto v0s = (TTree*)inputFile->Get(Form(
"%s/%s", dfName, v0Name.Data()));
206 if (v0s ==
nullptr) {
207 printf(
"%s table not found\n", v0Name.Data());
213 if (hasTrackQA && (trackQA = (TTree*)inputFile->Get(Form(
"%s/%s", dfName, trackQAName.Data()))) ==
nullptr) {
219 int trackIdxPos = 0, trackIdxNeg = 0;
220 std::vector<bool> keepV0s(trackExtraTree->GetEntries(),
false);
221 v0s->SetBranchAddress(
"fIndexTracks_Pos", &trackIdxPos);
222 v0s->SetBranchAddress(
"fIndexTracks_Neg", &trackIdxNeg);
223 auto nV0s = v0s->GetEntriesFast();
224 for (
int i{0};
i < nV0s; ++
i) {
226 keepV0s[trackIdxPos] =
true;
227 keepV0s[trackIdxNeg] =
true;
231 std::vector<bool> keepTrackQA;
233 keepTrackQA.assign(trackExtraTree->GetEntries(),
false);
234 trackQA->SetBranchAddress(
"fIndexTracks", &trackIdxPos);
235 for (
int i{0};
i < trackQA->GetEntries(); ++
i) {
236 trackQA->GetEntry(
i);
237 keepTrackQA[trackIdxPos] =
true;
241 std::vector<int> acceptedTracks(trackExtraTree->GetEntries(), -1);
242 std::vector<bool> hasCollision(trackExtraTree->GetEntries(),
false);
244 uint8_t tpcNClsFindable = 0;
245 bool bTPClsFindable =
false;
246 uint8_t ITSClusterMap = 0;
247 UInt_t ITSClusterSizes = 0;
248 bool bITSClusterMap =
false;
249 bool bITSClusterSizes =
false;
250 uint8_t TRDPattern = 0;
251 bool bTRDPattern =
false;
253 bool bTOFChi2 =
false;
256 TBranch* br =
nullptr;
257 TIter next(trackExtraTree->GetListOfBranches());
258 while ((br = (TBranch*)next())) {
259 TString brName = br->GetName();
260 if (brName ==
"fTPCNClsFindable") {
261 trackExtraTree->SetBranchAddress(
"fTPCNClsFindable", &tpcNClsFindable);
262 bTPClsFindable =
true;
263 }
else if (brName ==
"fITSClusterMap") {
264 trackExtraTree->SetBranchAddress(
"fITSClusterMap", &ITSClusterMap);
265 bITSClusterMap =
true;
266 }
else if (brName ==
"fITSClusterSizes") {
267 trackExtraTree->SetBranchAddress(
"fITSClusterSizes", &ITSClusterSizes);
268 bITSClusterSizes =
true;
269 }
else if (brName ==
"fTRDPattern") {
270 trackExtraTree->SetBranchAddress(
"fTRDPattern", &TRDPattern);
272 }
else if (brName ==
"fTOFChi2") {
273 trackExtraTree->SetBranchAddress(
"fTOFChi2", &TOFChi2);
280 if (!bTPClsFindable || !bTRDPattern || !bTOFChi2 ||
281 (!bITSClusterMap && !bITSClusterSizes)) {
282 printf(
" *** FATAL *** Branch detection failed in %s for trackextra.[(fITSClusterMap=%d,fITSClusterSizes=%d),fTPCNClsFindable=%d,fTRDPattern=%d,fTOFChi2=%d]\n", dfName, bITSClusterMap, bITSClusterSizes, bTPClsFindable, bTRDPattern, bTOFChi2);
287 int fIndexCollisions = 0;
288 track_iu->SetBranchAddress(
"fIndexCollisions", &fIndexCollisions);
291 auto entries = trackExtraTree->GetEntries();
293 for (
int i = 0;
i < entries;
i++) {
294 trackExtraTree->GetEntry(
i);
295 track_iu->GetEntry(
i);
298 hasCollision[
i] = (fIndexCollisions >= 0);
301 if (tpcNClsFindable > 0 && TRDPattern == 0 && TOFChi2 < -1. &&
302 (!bITSClusterMap || ITSClusterMap == 0) &&
303 (!bITSClusterSizes || ITSClusterSizes == 0) &&
304 (!hasTrackQA || !keepTrackQA[
i]) &&
312 for (
auto key2 : *treeList) {
313 TString treeName = ((TObjString*)key2)->GetString().Data();
317 if (outputDir ==
nullptr) {
318 outputDir = outputFile->mkdir(dfName);
319 printf(
"Writing to output folder %s\n", dfName);
323 auto inputTree = (TTree*)inputFile->Get(Form(
"%s/%s", dfName, treeName.Data()));
324 printf(
" Processing tree %s with %lld entries with total size %lld\n", treeName.Data(), inputTree->GetEntries(), inputTree->GetTotBytes());
325 auto outputTree = inputTree->CloneTree(0);
326 outputTree->SetAutoFlush(0);
328 std::vector<int*> indexList;
329 std::vector<char*> vlaPointers;
330 std::vector<int*> indexPointers;
331 TObjArray* branches = inputTree->GetListOfBranches();
332 for (
int i = 0;
i < branches->GetEntriesFast(); ++
i) {
333 TBranch* br = (TBranch*)branches->UncheckedAt(
i);
334 TString branchName(br->GetName());
335 TString tableName(
getTableName(branchName, treeName.Data()));
337 if (!tableName.EqualTo(
"O2track")) {
341 if (((TLeaf*)br->GetListOfLeaves()->First())->GetLeafCount() !=
nullptr) {
342 printf(
" *** FATAL ***: VLA detection is not supported\n");
344 }
else if (branchName.BeginsWith(
"fIndexSlice")) {
347 vlaPointers.push_back(
reinterpret_cast<char*
>(
buffer));
348 inputTree->SetBranchAddress(br->GetName(),
buffer);
349 outputTree->SetBranchAddress(br->GetName(),
buffer);
351 indexList.push_back(
buffer);
352 indexList.push_back(
buffer + 1);
353 }
else if (branchName.BeginsWith(
"fIndex") && !branchName.EndsWith(
"_size")) {
356 indexPointers.push_back(
buffer);
358 inputTree->SetBranchAddress(br->GetName(),
buffer);
359 outputTree->SetBranchAddress(br->GetName(),
buffer);
361 indexList.push_back(
buffer);
365 const bool processingTracked = treeName.BeginsWith(
"O2tracked");
366 const bool processingTrackQA = treeName.BeginsWith(
"O2trackqa");
367 const bool processingTracks = treeName.BeginsWith(
"O2track") && !processingTracked && !processingTrackQA;
368 const bool processingAmbiguousTracks = treeName.BeginsWith(
"O2ambiguoustrack");
370 auto entries = inputTree->GetEntries();
371 for (
int i = 0;
i < entries;
i++) {
372 inputTree->GetEntry(
i);
373 bool fillThisEntry =
true;
375 if (processingTracks) {
376 if (acceptedTracks[
i] < 0) {
377 fillThisEntry =
false;
381 for (
const auto& idx : indexList) {
382 int oldTrackIndex = *idx;
385 if (oldTrackIndex >= 0) {
386 if (acceptedTracks[oldTrackIndex] < 0) {
387 fillThisEntry =
false;
389 *idx = acceptedTracks[oldTrackIndex];
396 if (processingAmbiguousTracks) {
397 if (hasCollision[
i]) {
398 fillThisEntry =
false;
407 if (entries != outputTree->GetEntries()) {
408 printf(
" Reduced from %lld to %lld entries\n", entries, outputTree->GetEntries());
411 static const std::array<TString, 4> checkNames{
"O2track",
"O2trackextra",
"O2trackcov",
"O2ambiguoustrack"};
412 std::vector<bool> checks(checkNames.size(),
false);
413 for (
size_t i{0};
i < checkNames.size(); ++
i) {
414 if (tableName.EqualTo(checkNames[
i])) {
418 if (std::none_of(checks.begin(), checks.end(), [](
bool b) { return b; })) {
420 printf(
" -> Reduction is not expected for this tree!\n");
427 for (
auto&
buffer : indexPointers) {
430 for (
auto&
buffer : vlaPointers) {
451 printf(
"Removing incomplete output file %s.\n", outputFile->GetName());
452 gSystem->Unlink(outputFile->GetName());
459 auto sBefore = inputFile->GetSize();
460 auto sAfter = outputFile->GetSize();
461 if (sBefore <= 0 || sAfter <= 0) {
462 printf(
"Warning: Empty input or output file after thinning!\n");
465 auto spaceSaving = (1 - ((double)sAfter) / ((double)sBefore)) * 100;
466 printf(
"Stats: After=%lld / Before=%lld Bytes ---> Saving %.1f%% diskspace!\n", sAfter, sBefore, spaceSaving);
467 printf(
"Timing: CPU=%.2f (s); Real=%.2f (s)\n", clock.CpuTime(), clock.RealTime());
468 printf(
"End of AOD thinning.\n");