]> Creatis software - clitk.git/blob - tools/clitkMergeRootFiles.cxx
Merge branch 'master' into clitkImage2Dicom
[clitk.git] / tools / clitkMergeRootFiles.cxx
1 /**
2    =================================================
3    * @file   clitkMergeRootFiles.cxx
4    * @author David Sarrut <david.sarrut@creatis.insa-lyon.fr>
5    * @author Brent Huisman <brent.huisman@insa-lyon.fr>
6    * @date   06 May 2014
7    *
8    * @brief
9    *
10    =================================================*/
11
12 #include "clitkMergeRootFiles_ggo.h"
13 #include "clitkCommon.h"
14 #include "GateMergeManager.hh"
15 #include <string>
16 #include <TROOT.h>
17 #include <TPluginManager.h>
18 #include <TFile.h>
19 #include <TKey.h>
20 #include <TFileMerger.h>
21 #include <TTree.h>
22 #include <TChain.h>
23 #include <TH1.h>
24 #include <TH2.h>
25 #include <iostream>
26 #include <TApplication.h>
27
28 using std::endl;
29 using std::cout;
30
31 struct PetInputFile {
32     string filename;
33     double mean_time;
34 };
35
36 bool sort_pet_input_file(const PetInputFile &a, const PetInputFile &b) {
37     return a.mean_time < b.mean_time;
38 };
39
40 typedef std::vector<PetInputFile> PetInputFiles;
41
42 //-----------------------------------------------------------------------------
43 int main(int argc, char *argv[]) {
44         //this fixes a bug in TFileMerger. See http://root.cern.ch/phpBB3/viewtopic.php?t=18016.
45     TApplication app("app", 0, 0);
46
47     gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo", "*",
48                                           "TStreamerInfo", "RIO", "TStreamerInfo()");
49
50     // init command line
51     GGO(clitkMergeRootFiles, args_info);
52
53     // Check parameters
54     if (args_info.input_given < 2) {
55         FATAL("Error, please provide at least two inputs files");
56     }
57
58     /* The following block does some bookkeeping necesary for files originating from a pet simulation.
59     Seems fixing some timing info, for coincidences between files perhaps.
60     It seems the files are later on reopened and merged, if the conditions were met.
61     It's not required for merging other .root files.
62     GateMergeManager reportedly exists specifically for the purpose of merging pet simulations. */
63     {
64         // Detect Pet output
65         bool all_pet_output = true;
66         PetInputFiles pet_input_files;
67         for (uint i = 0; i < args_info.input_given; i++) {
68             const char *filename = args_info.input_arg[i];
69             PetInputFile input_file;
70             input_file.filename = filename;
71             TFile *handle = TFile::Open(filename, "READ");
72             assert(handle);
73             TTree *hits = dynamic_cast<TTree *>(handle->Get("Hits"));
74             TTree *singles = dynamic_cast<TTree *>(handle->Get("Singles"));
75             const bool is_pet_output = (hits != NULL) && (singles != NULL);
76             cout << "testing " << filename << " is_pet_output " << is_pet_output;
77
78                 //TTree *histos = dynamic_cast<TH1F *>(handle->Get("histo;1"));
79                 //const bool is_hist_output = (histos != NULL);
80
81
82             if (is_pet_output) {
83                 double time;
84                 double time_accum = 0;
85                 singles->SetBranchAddress("time", &time);
86                 size_t total = singles->GetEntries();
87                 for (size_t kk = 0; kk < total; kk++) {
88                     singles->GetEntry(kk);
89                     time_accum += time;
90                 }
91
92                 input_file.mean_time = time_accum / total;
93                 pet_input_files.push_back(input_file);
94                 cout << " mean_time " << input_file.mean_time;
95             }
96
97             cout << endl;
98
99             handle->Close();
100             delete handle;
101             //bitwise on booleans?
102             all_pet_output &= is_pet_output;
103         }
104         cout << "all_pet_output " << all_pet_output << endl;
105
106         if (all_pet_output) {
107             GateMergeManager manager(args_info.fastmerge_given, args_info.verbose_arg, true, 0, "");
108
109             cout << "sorting input file using singles time" << endl;
110             std::sort(pet_input_files.begin(), pet_input_files.end(), sort_pet_input_file);
111
112             Strings input_filenames;
113             for (PetInputFiles::const_iterator iter = pet_input_files.begin(); iter != pet_input_files.end(); iter++)
114                 input_filenames.push_back(iter->filename);
115
116             manager.StartMergingFromFilenames(input_filenames, args_info.output_arg);
117             //if the file was PET, then we're done.
118             return 0;
119         }
120     }
121
122     //if the file was not PET, but a generic Rootfile, use TFileMerger.
123
124     // Merge
125     TFileMerger *merger = new TFileMerger;
126     for (uint i = 0; i < args_info.input_given; i++) merger->AddFile(args_info.input_arg[i]);
127     merger->OutputFile(args_info.output_arg);
128     bool whatbool = merger->Merge();
129
130     cout << "whatbool: " << whatbool << " no more whatbool" << endl;
131     // this is the end my friend
132     return 0;
133     
134 }
135 //-----------------------------------------------------------------------------